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Chapter  1 
(J.  A.  Smith) 

EXECUTIVE  SUMMARY 

The  primary  objective  of  the  work  undertaken  in  this  project  was  to 
review  the  mathematical  structure  and  physical  concepts  inherent  in 
electromagnetic  terrain  feature  canopy  models.  One  goal  of  such  a 
review  was  to  determine  if  a  factoring  of  canopy  electromagnetic  inter¬ 
actions  into  separate  geometric  and  radiative  terms  could  be  achieved. 
A  second  goal  was  the  extension  of  modeling  concepts  to  multidimensional 
problems  either  spatially  (for  the  reflective  regime)  or  temporally  (for 
the  thermal  regime). 

We  have  been  successful  in  the  primary  objective  of  the  project  and 
in  the  analysis  of  the  factoring  of  geometric  dependent  terms  from 
individual  optical  properties.  In  addition  we  have  expanded  our  reper¬ 
toire  of  calculational  procedures  for  general  terrain  radiative  transfer 
calculations  and  have  addressed  several  issues  related  to  the  comparison 
of  models  and  experiments.  We  have  been  partially  successful  in  the 
analysis  of  multidimensional  problems. 

For  the  most  part,  the  major  results  of  this  study  have  been 
presented  in  the  reviewed  literature.  However,  Chapters  2-4,  discussed 
below,  present  a  detailed  discussion  of  the  mathematical  treatment  of 
canopy  biophysical  attributes  relative  to  the  calculation  of  canopy 
element  scattering  properties.  In  addition,  our  analysis  of  the  time- 
dependent  thermal  canopy  problem  has  not  yet  been  documented. 


Our  basic  modeling  approaches  in  both  the  thermal  and  reflective 
regime  are  given  in  Smith  (1980),  Appendix  Al,  Kimes  and  Smith  (1980), 


Appendix  A2,  Kimes,  Smith  and  Link,  (1981)  Appendix  A3  and  Cooper,  Smith 
and  Pitts  (1982),  Appendx  AS. 


The  factorization  of  the  thermal  model  into  geometrical  and 
optical-dependent  terms  is  given  in  Smith,  Ranson,  Nguyen,  Balick,  Link, 
Fritschen,  and  Hutchison  (1981),  Appendix  A4.  This  paper  also  indicates 
a  method  for  the  inverse  problem  relative  to  canopy  geometry.  The 
respective  factorization  for  the  reflective  problem  is  given  in  Cooper, 
Smith,  and  Pitts  (1982),  Appendix  AS. 

Several  issues  relative  to  the  comparison  of  models  with  field  data 
and  vice-versa  are  documented  in  Kimes,  Smith,  and  Ranson  (1980), 
Appendix  A6,  Kirchner,  Youkhana,  and  Smith  (1982),  Appendix  A7,  and 
Smith  (1982),  Appendix  A8. 

Finally,  in  collaborative  fashion  our  thermal  modeling  work  has 
been  used  in  the  comparison  of  various  modeling  abstractions  for  soil 
and  soil/vegetation  composites  in  Balick,  Scoggins,  and  Link  (1981), 
Appendix  A9,  and  in  three  dimensional  modeling  concepts  incorporated  in 
Kimes  and  Kirchner  (1982),  Appendix  AlO. 

The  remaining  chapters  of  the  report  all  deal  with  reflective 
modeling  of  plane  parallel  vegetation  canopies  with  particular  attention 
to  calculation  of  media  scattering  properties  in  terms  of  the  distribu¬ 
tion  of  canopy  scattering  facets.  The  general  problem  studied  is  the 
application  of  the  Duntley  equations  to  solution  of  the  vegetation 
radiative  transfer  problem.  Chapter  4  presents  a  rather  exhaustive 
analysis  of  the  problem.  Expressions  are  derived  for  all  five  of  the 
Duntley  equation  coefficients  in  terms  of  general  biophysical  attributes 
of  the  canopy  and,  in  particular,  the  complete  leaf  slope  distribution. 
A  general  solution  approach  is  outlined,  but  a  more  simplified  analysis 


is  presented  for  the  special  case  of  a  canopy  with  symmetrical  upper  and 
lower  leaf  optical  properties.  The  results  of  Chapter  4  are  of  theo¬ 
retical  interest  but  of  limited  practical  application  because  of  the 
general  complexity  of  the  analysis  and  lack  of  explicit  analytic  solu¬ 
tions.  We  have,  however,  used  some  of  the  expressions  in  more  tractible 
formulations,  e.g..  The  Adding  Method  of  Appendix  A5. 

Of  more  practical  consideration  is  a  simplified  analysis  of  the 
problem  presented  in  Chapter  2  where  only  the  Duntley  equation  coeffi¬ 
cients  dealing  with  the  direct  flux  components  are  analysed  in  detail. 
Tractible  expressions  are  formulated  and  evaluated  for  a  number  of 
standard  theoretical  canopy  geometries.  A  Suits  prime  model  incorporat¬ 
ing  the  extinction  coefficient  calculated  using  the  complete  leaf  slope 
distribution  is  formulated  and  compared  with  the  standard  Suits  model 
which  assumes  only  horizontal  and  vertical  canopy  projections  in 
Chapter  3.  Comparisons  with  both  theoretical  and  field  collected  data 
as  well  as  with  the  Monte  Carlo  SRVC  model  developed  under  previous  Army 
Research  Office  funding  are  presented.  The  new  Suits  prime  model  yields 
results  comparable  to  the  Monte  Carlo  solution  but  with  corresponding 
execution  times  reduced  by  a  factor  of  several  hundred. 

At  various  times  the  following  individuals  participated  in  the 
project: 

Dr.  James  A.  Smith,  Principal  Investigator 

Dr.  R.  Laven 

Mr.  S.  Youkhana,  Graduate  Research  Assistant.  Mr.  Youkhana 
received  the  degree  Doctor  of  Philosophy  for  work  related  to 
project  objectives. 

Mr.  H.  Randolph,  Graduate  Research  Assistant. 

Mr.  C.  Nikolopoulos,  Graduate  Research  Assistant.  Mr.  Nikolopoulos 
received  the  degree  Doctor  of  Philosophy  but  for  thesis 
research  unrelated  to  project  objectives. 
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A  list  of  all  publications  related  to  project  objectives  is  given 
in  Appendix  A  together  with  copies  of  journal  abstracts. 
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Chapter  2 

(S.  Youkhana  and  J.  A.  Smith) 
DIRECT  FLUX  EXTINCTION  COEFFICIENT 


Introduction 

Light  is  one  of  the  important  environmental  condition  parameters 
and  to  some  degree  a  critical  factor  to  the  green  plant  growth  and  its 
distribution.  Plant  growth  and  productivity  depend  on  the  quality  of 
solar  radiation  available  and  the  plant's  ability  to  utilize  visible 
light. 

Shading  caused  by  upperstory  (interception  of  light  by  dominant 
plants)  will  affect  the  understory  growth  (saplings,  seedlings,  and 
grasses).  Thus,  a  decrease  in  the  amount  of  light  for  understory  is 
mainly  due  to  absorption  of  light  by  leaves.  Absorption  of  light  is 
the  main  factor  governing  plant  photosynthesis  and  transpiration. 
Norman  (1980)  stated  that  vegetation  leaves  absorb  about  85  percent  of 
incident  visible  radiation  for  photosynthesis.  Understanding  the 
relation  between  photosynthesis  and  light  intensity  will  be  an  impor¬ 
tant  factor  to  understand  and  relate  dry  matter  to  other  community 
characteristics.  One  of  the  most  important  parameters  for  analysis  of 
production  rates  of  plant  communities  is  the  photosynthetic  rate 
(McCree  and  Loomis,  1969). 

As  mentioned,  interception  and  absorption  will  attenuate  light  as 
it  passes  through  the  canopy.  Many  researchers  experimentally  examined 


light  depletion  in  a  plant  community  as  a  function  of  plant  height  to 
examine  the  following: 

-  community  productivity 

-  morphology  of  understory 

*  species  that  dominate  as  succession  proceeds  toward 
climax 

-  invasion  and 

-  for  the  purpose  of  planning  and  vegetation  management. 

Community  productivity  was  the  main  purpose  of  investigators  since  it 
has  been  known  that  photosynthetic  rate  (CO2  uptake)  is  proportional  to 
the  amount  of  light  absorbed  up  to  a  compensation  point  (Donhoe,  1981). 
Some  effort  has  been  made  to  examine  leaf  morphology  at  different 
positions  in  plant  canopy.  Understanding  the  relationship  between 
light  intensity  and  community  characteristics,  a  proper  management 
decision  can  be  made.  With  proper  light  models,  one  can  make  decisions 
as  to  the  ideal  planting  pattern,  the  best  leaf  arrangement  and  the 
plant  properties  whi  should  be  changed  to  achieve  the  greatest  in¬ 
crease  in  water-use  efficiency  of  a  given  plant  species  (Lemeur  and 
Blad,  1974).  Thus,  all  field  data  previously  collected  was  by 
instantaneous  readings  inside  and  outside  plant  comnunities ,  and  then 
correlating  light  intensity  at  different  levels  to  photosynthetic  rate. 

Beginning  in  1950,  many  investigators  explained  light  attenuation 
at  different  canopy  heights  matheiaatically.  A  mathematical  model  is  an 
abstraction  of  the  real  world.  In  general,  for  the  model  to  be  appli¬ 
cable,  it  should  be  simple  to  understand  and  easy  to  handle.  However, 
in  order  that  the  mathematical  model  approximate  the  real  system 
behavior,  it  should  include  the  important  parameters.  This  means  that 
the  model  cannot  be  realistic  if  it  is  too  simple  and  its  adequacy  can 
only  be  judged  in  terms  of  the  purpose  for  which  it  has  been  created. 


Mathematical  models  for  the  extinction  coefficient  of  light 
through  a  plant  canopy  have  been  developed  over  a  period  of  time  by 
many  investigators  to  consider  both  diffuse  and  direct  light,  effects 
due  to  solar  zenith  angles,  and  canopy  characteristics.  One  of  the 
main  problems  in  modeling  the  extinction  coefficient  is  the  abstraction 
of  canopy  structure.  In  this  chapter,  the  extinction  coefficient  will 
be  formulated  by  projecting  canopy  elements  in  the  direction  of  light 
and  weighting  these  projections  by  the  leaf  distribution  function. 
This  approach  is  new  and  forms  the  basis  for  Suits  prime  (Suits') 
canopy  reflectance  model  discussed  later. 

Background 

Definitions  and  expressions  for  the  extinction  coefficient  of 
light  through  vegetation  canopies  have  been  changed  and  progressed  over 
a  period  of  time.  Projection  of  foliage  area  has  been  used  as  the 
basis  in  almost  all  mathematical  models  derived  for  the  extinction 
coefficient. 

Monsi  and  Saeki  (1953)  were  the  first  to  develop  theoretical  and 
experimental  methods  to  explain  and  understand  the  depletion  of  iso¬ 
tropic  diffuse  light  through  a  plant  community.  Their  experimental 
field  method  "stratifying  clip  technique"  related  observed  field  data 
of  canopy  characteristics  (accumulation  of  leaf  index  area)  to  light 
intensities  at  different  canopy  levels.  The  plotted  field  data  showed 
that  the  light  intensity  decreases  exponentially  as  leaf  index  area 
increases  above  a  certain  level.  Based  on  the  observed  relationship 
between  light  intensity  and  leaf  area  index  accumulation,  they 
formulated  their  expression  in  the  general  form  analogous  to 
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Beer* s-Lambert  law  of  monochromatic  flux  extinction  through  homogeneous 
absorbing  medium 

T  T  -KA 

I  =  I  e 
o 

where 

I  is  the  downward  intensity  of  isotropic  light  at  leaf  area 
index  A 

is  the  light  intensity  on  a  horizontal  surface  at  the 
canopy  top 

K  is  the  extinction  coefficient  of  light. 

Leaf  area  index  is  downward  accumulation  of  leaf  area 
(i.e.  A  =  A(x)dx 

t,  b,  X  are  canopy  top,  canopy  bottom  and  increment  in  leaf  area 
respectively.  The  coefficient  of  light  attenuation  is  given  by  the 
authors  as: 


A  A. 

K  =  ^  in  (1  -  ^) 
h  s 

where 

Aj^  is  the  mean  projection  of  actual  leaf  area  index  (A)  on 
horizontal  surface 

A^  is  the  soil  surface  area  above  which  foliage  area  is  measured. 
The  mean  projection  of  leaf  area  on  a  horizontal  surface  is  estimated 
as: 


A 


h 


A’sinB 

2  *  -1  2 
-A{sin6  sin  (tan6  cot6  )  cos6  tanO  (1  -  tan6 
71  *  a  a  r  a  r  a 


®r^®a 

cole  )^}  e  >e 

r"^  '  r  a 


They  claimed  that  there  was  an  agreement  between  measured  field  data 
and  model  prediction  of  extinction  coefficient. 
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Dewit  (1965)  extended  and  generalized  Monsi  and  Saeki's  theory  of 
light  attenuation  to  include  canopy  geometry.  His  definition  of  the 
extinction  coefficient  for  a  canopy  of  randomly  oriented  leaves  is  the 
ratio  of  foliage  unit  area  to  the  soil  unit  area  both  projected  in  the 
direction  of  light. 

Isobe  (1962)  used  cosine  direction  to  introduce  the  effect  of 
source  position  in  his  model  for  predicting  the  extinction  coefficient 
of  light.  He  also  proposed  a  distinct  model  for  leaves  and  for  stems 
in  the  same  canopy.  The  average  of  light  attenuation  through  canopy  is 
the  sum  of  the  contributions  from  both  leaves  and  stems.  As  in  Monsi 
and  Saeki's  work  (1953),  the  author's  model  will  predict  an  extinction 
coefficient  ranging  from  zero  to  unity  for  horizontally  oriented  leaves 
when  the  energy  source  is  at  horizon  and  at  nadir  respectively. 

-n'sine  0  se 

a  r  a 

^  ”  ■§{ (<l>”5)n' sin6  +  Jl'cosO  sin(|>}  6  >6 

^  d  Si  IT  3 

((i  =  cos  ^(■9t  tan  0  )  and  ?^(I)<71 

where  £'  and  n'  are  transformed  direction  cosines. 

Anderson  (1966)  reviewed  previous  work  on  light  attenuation 
through  plant  canopies  and  also  modified  the  definition  of  the  extinc¬ 
tion  coefficient  to  include  the  source  illumination  angle.  She  defined 
the  extinction  coefficient  in  terms  of  mean  foliage  area  projected  onto 
a  plane  normal  to  the  direction  making  an  angle  p  with  the  horizontal 
and  the  actual  leaf  area  index.  Both  Anderson  (1966)  and  Warren  Wilson 
(1967)  showed  the  extinction  coefficient  of  light  for  a  stand  of  uni¬ 
form  (constantly)  inclined  leaves  to  be  of  the  following  relation: 
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sin6^  ®r~®a 

^  ~  sin0  {l+§(tan(t>-(t») }  0  >0 

a  n  ^  ^  r  a 

(|>  =  cos  ^(tan0  cot0  ) 
a  r 

This  result  has  also  been  shown  by  Lupton  (1972).  In  1969,  Anderson 
derived  the  extinction  coefficient  for  randomly  oriented  vertical 
leaves  as: 

K  =  §tan0  0=0 

nr  a 

Cowan  (1968)  expressed  the  extinction  coefficient  as  the 
projection  mean  of  leaf  area  in  a  plane  perpendicular  to  incident 
radiation.  The  author  extended  Anderson’s  (1966)  idea  to  introduce 
three  different  idealized  abstract  plant  community  "stands"  as  follows: 

il.O  Canopy  with  horizontally  oriented  leaves 

^sec0^  Canopy  with  randomly  oriented  leaves 

2 

-tan0^  Canopy  with  vertical  leaves  randomly  oriented 

Recently,  investigators  have  collected  field  data  to  examine  the 
variation  in  extinction  coefficient  with  different  cover  types.  Miller 
(1967)  and  Landesberg  (1973)  collected  field  data  under  forest  type 
cover  to  examine  the  relation  between  the  extinction  coefficient  and 
source  elevation  and  mean  leaf  inclination  angles.  Similar  work  has 
been  done  for  agricultural  crops  (oats,  beans,  sunflower,  and  corn)  by 
Iinpens  and  Lemeur  (1969). 

Suits  (1972)  utilized  an  extinction  coefficient  as  a  function  of 
source  angle  0^  and  canopy  characteristics  in  his  reflectance  model. 
He  approximated  canopy  elements  by  projecting  them  orthogonally  in 
horizontal  and  vertical  planes. 
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K  =  H  +  §  V  tane 

71  r 

n 

H  =  A  f(e  )  cose  d0 
•'o  a  a  a 

n 

V  =  A  f(0  )  sine  de 

•'o  a  a  a 

H  and  V  are  orthogonal  projections  of  canopy  components  in 
horizontal  and  vertical  planes. 

In  this  chapter  a  more  complete  development  of  the  extinction 
coefficient  incorporating  the  complete  leaf  angle  distribution  is 
developed  and  contrasted  with  Suits  results  for  several  theoretical  and 
measured  canopies.  The  new  coefficient,  termed  K' ,  is  subsequently 
incorporated  into  the  Suits  model  to  yield  a  so-called,  Suits'  model. 
Model  comparisons  are  made  in  the  next  chapter. 

Approach 

The  direct  flux  extinction  coefficient  K  can  be  defined  as  the 
amount  of  power  (flux)  absorbed  (attenuated)  per  unit  path  length.  The 
amount  of  energy  absorbed  by  canopy  leaves  is  proportional  to  the 
amount  intercepted  by  those  elements.  The  amount  of  power  intercepted 
by  individual  canopy  elements  is  proportional  to  the  projection  of  the 
leaf  segment  in  the  direction  of  the  beam.  Projection  of  the  foliage 
area  in  the  direction  of  the  source  will  be  considered  as  the  basis  in 
this  approach.  Mathematical  equations  will  be  derived  to  predict  the 
direct  flux  extinction  coefficient  as  a  function  of  illumination  source 
and  leaf  distribution  function. 

Let  A  be  the  vector  area  of  an  individual  canopy  element  and  R  be 
the  unit  vector  in  the  direction  of  energy  source.  These  two  vectors 
can  be  written  in  cartesian  coordinates  as  following: 
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where 


A  =  Z  A.  e. 
i=l  ^  " 


R  =  Z  R.  8. 
1  1 

1 


A  •  A  =  A^ 

R  •  R  =  unity 


i  is  the  X,  Y,  and  Z  direction  respectively  and  is  a  unit  vector  in 
the  X,  Y,  and  Z  directions.  Both  vectors  A  and  R  can  be  written  in  the 
spherical  polar  (system)  coordinates  (Fig.  1)  where 
X  =  sin6  cos(ti 
Y  =  sin0  sin(t)  and 
Z  =  COS0 

A  =  A  8  +  A  8  +  A  8 
XX  y  y  z  2 


A  =  A  {sin0  cosib  e  +  sin0  sind)  e  +  cos0  e  } 
a  ^a  X  a  ^a  y  a  z 


(1) 


M  ^  ^  I 

R  =  {sin0  cos<t»  8  +  sin0  sin*  8  +  cos0  8  }  I 

'  r  * r  X  r  'r  y  r  f 

Where  0  and  (|i  are  zenith  and  azimuth  angles  respectively  of  each 

vector. 

The  power  intercepted  by  an  individual  canopy  element  A  from  the 
source  in  the  direction  R  is  proportional  to 


A  •  R  =  A  cos  (A  •  R) 

A  •  R  =  Z  A.R. 

1  1 

1 

A  •  R  =  A^R^  +  AyRy  +  A^R^  (2) 

Substituting  equation  1  into  2  will  yield 


i 

I 

i 


A  *  R  =  |a|  {sin6  cos(|)  sin6  cos^  +  sin6  sin<{)  sin6  sla(|> 


+  cose  cose  } 
a  r 


A  .  R  =  |A|  {sine^  sin6^  (cos<l)^  cos(|>^  +  sin<|)^  sin(|i^) 


+  cose  cose  } 
a  r' 


from  trigonometry  cos  ((|>  ~  (|>  )  =  cos$  cosiji  +  sin<l)  sincj) 

2k  T  dr  dr 


A  •  R  =  (a  I  {sine  sine  cos  (<)•  -(J>  )  +  cose  cose  }  (3) 

dr  dr  dr 

=  F((t)^  - 

For  simplicity,  let  (|»  =  (t»  -  (j»  which  is  the  azimuthal  angle  difference 

d  r 

between  source  direction  and  canopy  element. 

The  average  power  intercepted  per  canopy  element  A  (denoted 
<A  can  be  determined  from  the  total  interception  by  collection  of 

elements  randomly  distributed,  all  of  which  are  characterized  by  the 
same  element  zenith  angle  and  assuming  azimuthally  isotropic  leaves. 


=  /o  A-Rd4./;Jd(t) 


where  0  S  0  <  7t 


The  integral  will  be  broken  down  into  two  parts  based  on  a  critical 
azimuthal  angle  difference  between  energy  source  and  canopy  element. 
This  angle  determines  when  the  illumination  angle  is  greater  than  the 
foliage  angle.  When  azimuthal  angle  difference  0  is  smaller  than  the 
critical  angle  0^,  the  beam  will  strike  the  upper  side  of  the  element. 
The  beam  will  hit  the  lower  part  of  the  foliage  element  when  0  is 
greater  than  0^,  and  the  beam  will  be  parallel  to  the  foliage  element 
when  0=0. 


-It- 


<A-R>.  =  ^  A{/ 
(|>  n  ‘•'o 


A*Rd0  -  j"  A*Rd(t)} 


-  A  { T  (sin6  sin0  cosd)  +  cos6  cos0  )d<b 
n  ‘•'o  a  r  a  r  ^ 


■  Jj.(sin0  sin0  cosdi  +  cos0  cos0  )d6 
•^6  a  r  ^  a  r  ^ 


0  TT 

=  -  A  {sin0  sin0  sin*  +  ((|)  ~  x)  cos0  cos0  }  (4) 

71  a  r  o  2  a  r 

The  complication  arises  because  the  critical  azimuth  angle  (|»^  depends 

on  both  0^  and  0^  (foliage  and  source  zenith  angles).  can  be  found 

by  investigating  the  sign  of  F(<J)  -•!►),  assuming  that  both  0  and  0 

3  1C  d  IT 

are  in  the  range  between  zero  and  ?  (0  ^  0  or  0  ^  as  following: 

a.  If  ^^0+0<7i 

2  a  r 

cos(0  +  0  )  S  0 
a  r 

F(())  -((»)  =  sin0  sin0  cos0  +  cos0  cos0  ^  0 

^3  ^  I*  3  JT  3  T 

coscf)  ^  -COt0  COt0 
o  a  r 

(b  =  cos  ^(-cot0  COt0  I 
o  a  r 

b.  If  O<0+0s5 

a  r  2 

cos(0  +0)^0 
a  r 

F((}i  -  (|(  )  =  sin0  sin0  cos(^  +  cos0  cos0  ^  0 

a  ^r  a  r  o  a  r 

cos<b  =  -cot0  cot0  %  -1 
o  a  r 

COS0  =  “1  and  (1)  =  7i 

For  <I>  >  <t»^  =>  F((^)  %  0  and  for  (f  ^  F((»  ^  0 

One  of  the  main  problems  in  this  formulation  is  how  to  introduce 
canopy  abstraction  and  how  to  project  the  foliage  area  index.  Suits 
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(1972)  approximated  component  leaf  areas  by  the  sum  of  orthogonally 
projected  areas  in  horizontal  and  vertical  planes.  More  generally, 
however,  normalized  leaf  distribution  functions  (Fig.  2)  should  be 
introduced  to  predict  the  average  power  intercepted  by  the  collection 
of  leaves  over  their  inclination  angles.  The  proportion  of  leaves 
having  a  slope  between  6  and  0  dS  is  given  by  f(6  )d0  then 

-  -  0  - 

<A*R>.  „  =  (-)^  A  ^  <A*R>x  f(0  )d8  +  si  <A-R>.  f(6  )d0  } 

(I),6  ‘•'o  ^  ^  a  a  ■^n  »  <(>  '  a  a^ 

a  2  ■  ‘’r 


5  -  fl 
0  0  0 

=  (-)^  A  {/  ^sin0  sin0  sin())  f(6  )d0 
n  •'o  a  r  ^o  a  a 


-  -  0  - 

+  *^((1)  -  ?)  COS0  COS0  f(6  )d0  +  si  sin6 

•'o  ^o  2  a  r  '  a  a  •'n  -  a 

2  “  ®r 


n 

sin0  sin*  f(0  )d0  +  sl  (4»„  "  ?)  cos0  cos6  f(0  )d6  } 

r  033/1^0^  d  L  a  a 

2  ■  ®r 


The  critical  azimuthal  difference  angle  is  equal  to  n  in  the  range 

where  foliage  inclinations  are  between  zero  and  ^  This  means  the 

first  term  in  the  above  relation  will  approach  zero  and  the  last  three 
terms  will  predict  the  average  power  intercepted  by  a  single  element. 

-  -  0  - 

<A*R>.  a  “  §  A  COS0  S^  *^cos0  f(0  )d0  +  (§)^  A  sin0  sl  sin0 

(b,0  n  r  •^o  a  a  a  n  r  ■'n  a 

2  ■  ®r 


sin<|t^  f(0a)d03  +  (^)^  A  cos0^  ((^^  -  cos0^  fCe^)d0^  (5) 

2  ■  ®r 


Cumulative  Density  Function 
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Fig.  2.  Leaf  Distribution  Functions. 
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The  extinction  coefficient  is  given  by: 

_  Average  power  intercepted  (absorbed) 
unit  path  length 

The  cosine  of  the  source  zenith  angle  will  correct  for  increased  path 
lengths  at  varying  sun  elevations  thus  the  general  form  is  given  by: 

-  -  6  5 

2  ■  ®r 
n 

f(e,)d0^  t  (2)2  A  ;2  (♦__  -  S)  cose,  £(6,)de,  (6) 

2  ■  ®r 

Results 

Normalized  leaf-slope  distribution  functions  should  be  introduced 
in  equation  (6)  to  evaluate  direct  flux  extinction  coefficient.  The 
following  theoretical  abstractions  have  been  used  to  calculate  attenua¬ 
tion  coefficient: 

Planophile  -  horizontal  leaves  are  most  frequent 
f(6,)  =  2  (1  ♦  cos  26,) 

Erectophile  -  vertical  leaves  are  most  frequent 
f(ej  =  I  (1  -  cos  20  ) 

o  rt  a 

Plagiophile  -  oblique  leaves  are  most  frequent 

=  I  (1  - 

O  O 

Extremophile  -  oblique  leaves  are  less  frequent 

f(6a)  ■  n  ^  ^®a^ 

Uniform  -  same  proportion  leaves  at  any  angle 

f(0_)  =  I 


Spherical  -  Relative  frequency  of  leaf  inclination  is  the  same 
as  relative  frequency  of  inclination  of  surface 
element  of  a  sphere. 

f(6  )  =  sin6 
a  a 

Measured  leaf-slope  distributions  for  a  variety  of  canopies  were 
also  used. 

Direct  flux  extinction  coefficient  derived  by  Suit  (1972)  using 
orthogonal  projection  approximation  is  of  the  following  relation 

K  =  {<A  •R>.  a  •R>.  «  }/cos6 

V  h  (fr.Qg  r 

n  n 

K  =  LAI{/^  cose  f(e  )d0  +  I  tan6  sine  f(e  )de  }  (7) 

and  the  expression  for  attenuation  coefficeint  using  projection  of 
elements  in  the  direction  of  beam  approach  is  as  following: 


-  -  0  5.0 

K'  =  lai{2  -  !)  cose,  f(9,)de,  ♦  I  t„e,  '  sine. 


sine,  f(e,)d9,  *l!l  (♦,  -  5)  cose,  f(6,)de, 

2  ■  ®r 


n 

+  §  tane  sine  sin<^  f(e  )de  } 

n  r  •'n  _  a  ^o  a  a’ 

2  ■  ®r 


(8) 


The  first  and  third  terms  in  equation  (8)  account  for  the  first  term  in 
equation  (7)  which  is  leaf  area  projected  in  the  horizontal  plane  X-Y 
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coordinate.  The  second  and  fourth  terms  in  equation  (8)  account  for 
the  second  term  in  equation  (7)  which  is  the  sum  of  leaf  area  projected 
in  vertical  cartesian  coordinates  Y-2  and  X-Z  planes. 

The  critical  azimuthal  angle  difference  between  source  and  leaf 

plane  ()>^  is  equal  to  the  7t  when  leaf  zenith  angle  is  in  the  range 

between  zero  and  ^-6  (0<6  This  implies  that  part  of 

Z  r  a  z  r 

the  contribution  due  to  vertical  components  (second  term  in  equation 
(8)  will  approach  zero,  while  this  contribution  exists  in  orthogonal 
projection.  In  the  region  where  0  >  ^  -  6  ,  the  critical  azimuthal 

angle  difference  will  be  in  the  range  ^  and  n  (^  <  <  n) .  In  this 

range,  sin(t>^  will  take  values  less  than  unity  and  the  contribution  will 
be  minimized  proportional  to  the  value  of  sin(|i^.  It  should  be  noted 
that  this  reduction  does  not  exist  in  orthogonal  method.  Over  all,  the 
orthogonal  projection  technique  overestimates  the  average  amount  of 
power  received  by  canopy  elements  compared  to  what  will  actually  be 
intercepted. 

To  calculate  specific  extinction  coefficients  for  different  canopy 
geometries,  equation  (8)  must  be  solved.  Integration  of  equation  (8) 
over  leaf  inclination  angles  is  quite  difficult  because  of  the  depen¬ 
dence  of  (b  on  both  9  and  0  respectively,  6  =  cos  ^(-cot0  cot0  ). 
To  perform  the  integration  and  to  obtain  fairly  accurate  results, 
numerical  solutions  are  needed.  Simpson's  composite  algorithm  numeri¬ 
cal  method  (Burden  et.  al.,  1978)  over  2m  subinterval  of  [a,b]  was 
used  to  solve  the  problem  and  evaluate  I  =  f(x)dx  where 

f(x)dx  =  |{f(a)  +  2  V  f(x,,)  +  4  Z  f(x,.  ,)  +  f(b)} 

a  J  j=l  1 
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Extinction  coefficients  have  been  calculated  using  both  approaches 
for  different  theoretical  vegetation  architectures  having  a  leaf  area 
index  of  unity  (Table  1)  and  for  different  field  crops  (wheat/May, 
wheat/ June,  blue  grama,  and  sorghum)  having  different  leaf  area  indices 
(Table  2).  The  percentage  differences  in  predicted  K  values  between 
both  approaches  have  been  plotted  versus  source  zenith  angles  for  the 
four  theoretical  canopy  geometries  (Fig.  3)  and  for  measured  blue  grama 
grass,  wheat  and  sorghum  shown  in  Fig.  4. 

From  Figures  3  and  4,  it  is  clear  that  the  percentage  differences 
in  K  is  a  function  of  two  factors: 
a.  Canopy  abstraction 

In  general  the  percentage  difference  in  K  between  the  two 
methods  is  a  minimum  for  extremophile  canopy  because  most  of 
canopy  leaves  are  not  inclined  but  are  at  both  extreme  ends 
(vertical  and  horizontal).  Thus,  the  orthogonal  projection 
method  can  approximate  this  canopy  fairly  well  due  to  the 
presence  of  a  higher  proportion  of  leaves  in  vertical  and 
horizontally  oriented  planes.  Bunnik  (1978)  stated  that  the 
Suits  canopy  model  is  an  extreme  example  of  an  extremophile 
canopy.  On  the  other  hand,  the  greatest  percentage  differ¬ 
ences  in  K  is  found  for  plagiophile  canopies  because  inclined 
leaves  in  this  abstraction  are  more  frequent. 

For  a  spherical  abstraction  which  approximates  most 
grasses  and  crops  (e.g.  blue  grama.  Smith  and  Oliver  1972; 
soybeans,  Ranson  1981;  sorghum,  agronomy  field  data  1981), 
the  percentage  difference  in  K  values  ranges  from  7  to  43 
with  an  average  of  25  perent.  The  difference  is  mainly  due 
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Table  1.  Suits  and  Suits'  model  predictions  of  K  and  percentage 

difference  between  model  predictions  for  four  theoretical 
canopy  abstractions. 


0 

s 

Planophile 

Extremophile 

Plagiophile 

Spherical 

K' 

5 

0.848 

0.597 

0.769 

0.502 

15 

0.849 

0.616 

0.680 

0.518 

25 

0.850 

0.655 

0.685 

0.552 

35 

0.855 

0.715 

0.703 

0.610 

45 

0.870 

0.799 

0.749 

0.707 

55 

0.908 

0.926 

0.852 

0.872 

65 

1.008 

1.146 

1.086 

1.183 

75 

1.317 

1.668 

1.703 

1.932 

85 

3.229 

4.448 

4.970 

5.737 

K 

5 

0.872 

0.627 

0.717 

0.544 

15 

0.921 

0.696 

0.795 

0.634 

25 

0.975 

0.771 

0.881 

0.733 

35 

1.038 

0.859 

0.982 

0.850 

45 

1.119 

0.972 

1.111 

1.000 

55 

1.235 

1.134 

1.296 

1.214 

65 

1.428 

1.405 

1.606 

1.572 

75 

1.857 

2.006 

2.292 

2.366 

85 

3.937 

4.918 

5.620 

6.215 

|(K’  -  K)|/K' 

X  100 

5 

2.8 

5.0 

6.8 

8.4 

15 

8.5 

13.0 

16.9 

22.4 

25 

14.7 

17.7 

28.6 

32.8 

35 

21.4 

20.1 

39.7 

39.3 

45 

28.6 

21.7 

48.3 

41.4 

55 

36.0 

22.5 

52.1 

39.2 

65 

41.7 

22.6 

47.9 

32.9 

75 

41.0 

20.3 

34.6 

22.5 

85 

21.9 

10.6 

13.1 

8.3 
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Table  2.  Suits  and  Suits'  model  predictions  of  K  and  percentage 

difference  between  model  predictions  for  three  different 
crops . 


0 

s 

May/Wheat 

June/Wheat 

B.  Grama 

Sorghum 

5 

0.7346 

K’ 

0.4836 

3.6595 

4.259 

15 

0.7482 

0.4900 

3.6993 

4.316 

25 

0.7727 

0.5091 

3.8386 

4.507 

35 

0.8179 

0.5417 

4.1249 

4.941 

A5 

0.8993 

0.5965 

4.6392 

5.766 

55 

1.0433 

0.6911 

5.5791 

7.218 

65 

1.3270 

0.8320 

7.4021 

10.017 

75 

2.0446 

1.3260 

11.8477 

16.610 

85 

5.7912 

3.7142 

34.7030 

49.822 

K 

5 

0.8003 

0.5283 

3.9694 

4.653 

15 

0.8937 

0.5887 

4.5182 

5.448 

25 

0.9963 

0.6551 

5.1214 

6.321 

35 

1.1174 

0.7335 

5.8328 

7.352 

45 

1.2726 

0.8339 

6.7445 

8.672 

55 

1.4942 

0.9772 

8.0466 

10.557 

65 

1.8649 

1.2171 

10.2252 

13.712 

75 

2.6866 

1 . 7487 

15.0532 

20.703 

85 

6.6709 

4.3265 

38.4643 

54.604 

|(K'  -  K)|/K'  X  100 

5 

8.9 

9.2 

8.5 

9.3 

15 

19.4 

20.1 

22.1 

26.2 

25 

28.9 

28.7 

33.4 

40.2 

35 

36.6 

35.4 

41.4 

48.8 

45 

41.5 

39.8 

45.4 

50.4 

55 

43.2 

41.4 

44.2 

46.3 

65 

40.5 

39.4 

38.1 

36.9 

75 

31.4 

31.9 

27.1 

24.6 

85 

15.2 

16.5 

10.8 

9.6 

LAI  for  May/Wheat  =  1.28 
LAI  for  June/Wheat  =  0.82 
LAI  for  3.  Grama  =  6.51 
LAI  for  Sorghum  =  8.50 
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to  an  overestimate  in  vertical  projections.  Verhoef  and 
Bunnik  (1981)  stated  that  Suits  model  overpredicts  K  mainly 
due  to  drastic  simplifications  of  vertical  projection  leaf 
area. 

b.  Sun  zenith  angles 

The  percentage  differences  are  minimum  at  low  and  high  source 
angles  as  shown  in  Figs.  3  and  4.  This  is  due  to  contribu¬ 
tions  of  projection  canopy  elements  into  horizontal  plane  at 
low  energy  source  zenith  angles  and  projection  into  vertical 
planes  at  higher  zenith  angles.  Also,  the  difference  is  less 
when  the  source  is  at  nadir  than  at  the  horizon.  This  result 
is  expected  because  of  the  overestimate  arising  from  vertical 
components  by  the  orthogonal  projection  method  at  low  sun 
elevations . 

Conclusions 

A  mathematical  relationship  based  on  projection  of  a  single 
foliage  element  in  the  direction  of  an  energy  source  has  been  derived 
to  calculate  the  direct  flux  extinction  coefficient.  Calculation  of 
light  attenuation  through  the  canopy  has  been  performed  for  different 
standard  vegetation  canopy  architectures  (planophile,  extremophile, 
plagiophile,  and  spherical)  and  for  blue  gramma  grass,  wheat  and  sor¬ 
ghum,  using  both  the  Suits  orthogonal  projection  method  and  the  projec¬ 
tion  in  the  direction  of  source  methods. 

In  general,  the  orthogonal  projection  method  overestimates  the 
light  attenuation  coefficient  for  all  canopy  abstractions.  The  extra 
contribution  is  mainly  due  to  an  overestimation  of  foliage  elements 
projected  in  the  vertical  planes.  The  largest  differences  between  the 


Relative  Difference  (%) 
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Fig.  3.  Percentage  Difference  in  K  given  by 
(K-K')/K'  X  100. 


Sun  Zenith  Angle  (deg) 


Fig.  4.  Percentage  Difference  in  K  given  by 
(K-K')/K'  X  100. 
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two  approaches  was  found  to  be  for  the  canopy  whose  leaves  are  mostly 
inclined  (plagiophile)  and  the  smallest  differences  was  found  to  be  for 
the  canopy  in  which  most  of  its  leaves  are  not  inclined,  but  are  at 
extreme  angles  (extremophile) . 

The  extinction  coefficient  of  light  through  vegetation  canopy  will 
determine  how  the  beam  will  be  depleted  with  depth.  The  intensity  of 
light  reaching  the  soil  surface  under  the  canopy  will  be  inversely 
proportional  to  the  light  extinction  coefficient.  This  means,  the 
contribution  to  overall  canopy  reflectance  from  the  soil  surface  under 
the  canopy  will  be  underestimated  in  the  orthogonal  projection  method 
because  the  light  will  fall  off  in  the  canopy  faster  as  compared  to  the 
leaf  angle  distribution  method.  As  will  be  seen  later,  when  the  K' 
extinction  coefficient  is  inserted  into  the  Suits  radiative  transfer 
equation,  agreement  of  this  model  with  measurements  is  improved. 
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Chapter  3 

(S.  Youkhana  and  J.  A.  Smith) 
MODEL  COMPARISONS 


Three  models  were  analyzed  in  this  study  including  the  original 
Suits  model  (Suits,  1972),  the  modified  Suits  prime  developed  as  part 
oi  this  study  (Chapter  2)  and  the  SRVC  model  developed  by  Smith  and 
Oliver  (1972).  Complete  descriptions  of  the  original  Suits  and  SRVC 
models  are  given  in  the  original  literature  citations.  Briefly,  the 
models  differ  as  follows: 

The  Suits  model  is  an  extension  of  the  Duntly  equation  approach  in 
which  the  coefficients  of  the  differential  equations  are  expressed  in 
terms  of  the  optical  and  geometric  properties  of  the  canopy.  The 
optical  properties  are  taken  to  be  the  soil  reflectance  and  the  leaf 
reflectance  and  transmittance  as  a  function  of  wavelength.  The 
geometry  of  the  canopy  is  described  in  terms  of  its  leaf  area  index  and 
the  average  horizontal  and  vertical  leaf  area  projections.  A  direct 
solar  source  is  assumed  with  the  capture  of  solar  energy  governed  by 
the  direct  solar  extinction  coefficient.  The  model  predicts  the  flux 
reflected  from  the  canopy  as  a  function  of  sensor  view  angle. 

The  Suits'  model  developed  in  this  study  is  virtually  identical 
with  the  Suits  model  except  in  one  very  important  detail;  namely,  the 
calculation  of  the  direct  flux  extinction  coefficient  as  derived  in 
Chapter  2.  The  key  canopy  abstraction  which  is  modified  for  the  Suits' 
model  is  the  consideration  of  the  complete  leaf  angle  distribution 
rather  than  only  the  projected  horizontal  and  vertical  components. 
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The  SRVC  model  uses  a  stochastic  approach  to  simulate  the  light 
interactions  within  a  canopy  and  the  subsequent  reflection  from  it. 
Identical  canopy  optical  properties  are  used  as  in  the  previous  two 
models,  i.e.,  soil  reflectance  and  leaf  reflectance  and  transmittance; 
however,  the  complete  leaf  angle  distribution  is  used.  Since  leaves  in 
a  vegetation  canopy  do  not  form  a  systematic  pattern,  the  model  permits 
the  utilization  of  user  specified  statistical  rules.  The  model  differs 
from  the  Suits  model  then  in  two  respects:  the  use  of  a  leaf  angle 
distribution  and  it  is  not  an  analytical  solution.  It  differs  from  the 
Suits'  model  only  by  its  Monte  Carlo  nature  rather  than  its  conceptual 
abstraction.  Consequently,  it  was  hypothesized  that  the  modified  Suits 
model,  the  Suits'  model,  would  agree  more  closely  with  the  SRVC  model. 

To  evaluate  this  hypothesis  several  sets  of  measured  field  data 
were  used  in  which  both  the  model  input  parameters  and  the  model  pre¬ 
dicted  values  were  estimated.  The  complete  set  of  input  data  are  given 
in  Tables  3  and  4,  and  used  to  drive  all  the  model  siumlations  used  in 
this  study.  The  wheat  data  were  collected  by  a  Colorado  State 
University/Texas  A  &  M  University  experimental  team  during  April  17, 
May  16,  and  June  13,  1976  as  part  of  the  National  Aeronautics  and  Space 
Administration  Large  Area  Crop  Inventory  Exp  jriment  (LACIE)  and 
extracted  from  Smith,  Berry,  and  Heimes  (1976).  The  study  site  was  in 
Finney  County,  Kansas. 

Similarly,  the  blue  grama  data  set  was  collected  during  August, 
1972  by  Colorado  State  University  personnel  at  the  Pawnee  National 
Grasslands  at  the  Intensive  Study  Site  of  the  International  Biological 
Program  and  were  extracted  from  Smith  and  Oliver  (1972).  An  additional 
advantage  of  using  both  this  and  the  previous  data  set  was  that  the 
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Table  3.  Vegetation  canopy  characteristics  (leaf  area  index  and 

orthogonal  projection),  leaf  optical  properties  and  soil 
reflectance. 


Vegetation 

LAI 

H 

V 

\ 

P 

T 

Ps 

Sorghum 

8.50 

4.27 

6.92 

660 

0.108 

0.025 

0.085 

830 

0.405 

0.483 

0.263 

Wheat/ April 

0.87 

0.51 

0.59 

550 

0.078 

0.078 

0.156 

650 

0.067 

0.067 

0.218 

750 

0.331 

0.331 

0.256 

950 

0.480 

0.480 

0.311 

Wheat/May 

1.23 

0.76 

0.81 

550 

0.076 

0.076 

0.156 

650 

0.045 

0.045 

0.218 

750 

0.381 

0.381 

0.256 

950 

0.450 

0.450 

0.311 

Wheat/ June 

0.82 

0.50 

0.53 

550 

0.168 

0.168 

0.156 

650 

0.206 

0.206 

0.218 

750 

0.437 

0.437 

0.256 

950 

0.480 

0.480 

0.311 

Blue  grama 

6.5 

3.71 

4.78 

400 

0.038 

0.014 

0.024 

450 

0.047 

0.022 

0.027 

500 

0.059 

0.034 

0.032 

550 

0.098 

0.095 

0.042 

600 

0.074 

0.064 

0.046 

650 

0.053 

0.035 

0.050 

700 

0.116 

0.083 

0.065 

750 

0.501 

0.277 

0.098 
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Table  4.  Probability  distribution  for  leaves. 


p 

April 

Wheat 

May 

June 

Blue  grama 

Sorghum 

0 

0.044 

0.031 

0.050 

0.0141 

0.0031 

5 

0.044 

0.029 

0.050 

0.0197 

0.0049 

10 

0.044 

0.046 

0.057 

0.0248 

0.0061 

15 

0.044 

0.058 

0.055 

0.0280 

0.0071 

20 

0.045 

0.056 

0.047 

0.0312 

0.0940 

25 

0.046 

0.054 

0.047 

0.0349 

0.0205 

30 

0.047 

0.053 

0.047 

0.0449 

0.0304 

35 

0.048 

0.052 

0.047 

0.0541 

0.0434 

40 

0.049 

0.056 

0.052 

0.0664 

0.0537 

45 

0.051 

0.060 

0.053 

0.0768 

0.0729 

50 

0.052 

0.059 

0.048 

0.0727 

0.0814 

55 

0.054 

0.059 

0.050 

0.0699 

0.0997 

60 

0.055 

0.063 

0.052 

0.0845 

0.1174 

65 

0.057 

0.062 

0.054 

0.0964 

0.1287 

70 

0.059 

0.059 

0.055 

0.0880 

0.1191 

75 

0.062 

0.055 

0.051 

0.0827 

0.0836 

80 

0.064 

0.051 

0.060 

0.0593 

0.0680 

85 

0.067 

0.050 

0.064 

0.0391 

0.0288 

90 

0.070 

0.046 

0.060 

0.0123 

0.0218 

Total 

1.00 

1.000 

1.000 

1.000 

1.000 

7 
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SRVC  model,  which  is  expensive  to  run,  had  already  been  applied  to 
these  canopies. 

Finally,  a  third  data  set  was  collected  by  the  author  over  a 
s<n  ghum  field  on  August  22,  1981  at  the  Colorado  State  University 
Agronomy  Farm.  The  Goddard  Space  Flight  Center  Mark-II  three  band 
hand-held  radiometer  was  used  to  obtain  canopy  reflectance  measure¬ 
ments.  An  ISCO  radiometer  was  used  to  estimate  leaf  optical  proper¬ 
ties.  The  canopy  was  sampled  to  estimate  leaf  area  index  using  a 
LICOR-3000  Optical  portable  area  meter.  While  care  was  taken  to  obtain 
a  representative  measurement  for  the  homogeneous  canopy,  an  accurate 
estimate  of  sampling  error  could  not  be  obtained  since  destructive 
sampling  was  prohibited.  Subsequent  sensitivity  analysis,  indicated 
that  this  input  parameter  was  not  critical  at  the  visible  wavelengths 
but  was  of  concern  at  infrared  wavelengths.  Consequently,  the  infrared 
data  were  dropped  from  further  consideration. 

Table  5  summarizes  all  of  the  model  predictions  for  each  of  the 
three  data  sets.  9^  is  the  instrument  viewing  direction  and  0^  the 
sun  angle  at  the  time  of  measurement.  The  model  predictions  are 
plotted  against  field  measurements  for  blue  grama  in  Figures  5  and  6 
for  two  different  view  angles.  Figures  7  and  8  summarize  the  data  for 
the  wheat  canopy  for  a  vertical  view  angle,  and  finally.  Figures  9  and 
10  show  the  behavior  of  the  models  for  two  different  sun  angles  and 
varying  view  angles  for  the  sorghum  canopy  at  the  visible  wavelength. 

As  a  general  observation  it  is  evident  that  the  Suits'  predictions 
are  consistently  more  in  agreement  with  the  measured  data  and  the  SRVC 
model  predictions  than  hypothesized.  The  Suits  model  systematically 
underestimates  canopy  reflectance  as  compared  to  the  measured  data  and 
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Bidirectional  Reflectance  Factor 
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Woveiength 

Fig.  5.  SRVC,  Suits,  and  Suits'  Model  Predictions  Compared  with  the 

Measured  Blue  Gramma  Canopy  Reflectance  for  Sensor  Polar 
Angle  of  5*  and  Sun  Zenith  Angle  of  35”. 


Bidirectional  Reflectance  Factor 


Wavelength 


Fig.  7.  SRVC,  Suits,  ond  Suits'  Model  Predictions 
Compared  with  the  Measured  Wheat  Canopy 
Reflectance  for  Sensor  Polar  Angle  of  5* 
and  Sun  Zenith  Angle  of  40.5**. 


550  650  750  950 

Wavelength 

Fig.  8.  SRVC,  Suits,  and  Suits'  Model  Predictions 
Compared  with  the  Measured  Wheat  Canopy 
Reflectance  for  Sensor  Polar  Angle  of  5® 
and  Sun  Zenith  Angle  of  38.1”. 
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this  in  consistent  with  the  theoretical  calculations  presented  in 
Tables  1  and  2  of  Chapter  2  hich  indicate  that  this  model  over¬ 
estimates  the  direct  flux  extinction  coefficient.  Consequently,  less 
flux  will  be  reflected  from  the  underlying  soil  layer.  Less  flux  will 
also  be  reflected  from  the  leaves  since  absorption  is  assumed  greater. 

More  specifically,  it  is  seen  from  Table  5  that  the  difference 
between  the  Suits  and  Suits'  model  for  blue  for  blue  grama  is  between 
25  and  40  percent;  for  wheat  the  differences  are  on  the  order  of  20 
percent  and  about  30  percent  for  sorghum.  Generally,  the  results  are 
insensitive  to  soil  reflectance.  Leaf  optical  properties  should  be 
measured  to  the  same  degree  of  precision  as  the  overall  canopy 
reflectance  and  leaf  area  index  should  be  determined  to  within 
approximately  50  percent  of  its  nominal  value.  The  most  stringent 
requirement  experimentally  is  the  estimation  of  leaf  area  index  for  low 
plant  cover  situations.  This  could  be  a  limiting  factor,  therefore, 
for  the  wheat  data  set. 

As  a  further  general  comparison  test  of  the  three  models  a  simple 
Chi  Square  test  between  each  of  the  model  predictions  and  the  measured 
field  data  was  performed.  The  value  for  the  SRVC  model  was  23.6;  for 
the  Suits'  model,  28.8,  and  for  the  original  Suits  model,  40.3.  The 
corresponding  level  of  significance  with  29  degrees  of  freedom  for 
these  values  correspond  to  75  percent,  50  percent  and  10  percent 
respectively  in  general  agreement  with  the  observations  above. 


Chapter  4 

(C.  Nikolopoulos  and  J.  A.  Smith) 

A  MORE  COMPLETE  ANALYSIS  OF  THE  DUNTLEY  EQUATIONS  AND  THEIR 
APPLICATION  FOR  DESCRIBING  RADIATIVE  TRANSFER  IN 
VEGETATION  TERRAIN  ELEMENTS 

Statement  of  the  Problem  and  Plan  of  Attack 

The  Duntley  differential  equations  for  direct  and  diffuse  flux 
extinction  offer  a  tractible,  analytic  solution  approach,  for  radiative 
transfer  in  plane  parallel  media  whose  optical  properties  can  be  simply 
characterized.  Diffusing  plates  are  a  classic  example.  Vegetation 
canopies,  however,  are  more  complex  in  that  the  individual  scattering 

elements;  i.e,,  leaves,  asymmetrically  scatter  radiation  both  as  a 

function  of  their  intrinsic  optical  reflectance  and  transmittance  and 

their  orientations.  The  original  Suits  model  discussed  in  Chapter  2 
applied  the  Duntley  equations  to  vegetation  media  but  in  order  to 
develop  tractible  expressions  for  the  Duntley  coefficients,  abstracted 
the  canopy  elements  as  possessing  only  vertical  and  horizontal  projec¬ 
tions,  essentially  a  Dirac  delta  function  approximation  of  the  leaf 
slope  distribution.  The  Suits  prime  (Suits’)  model  developed  in 
Chapter  2  extended  the  calculations  to  include  consideration  of  the 

complete  leaf  slope  distribution  in  the  derivation  of  the  direct  flux 
extinction  coefficient  and  projected  view  angle.  As  was  seen  in 
Chapter  3,  this  modification  made  a  considerable  change  in  original 
Suits  model  predictions,  generally  an  improvement  by  a  factor  between  20 
and  40  percent. 

In  this  chapter,  the  complete  analysis  of  the  Duntley  equation 
coefficients  is  carried  out  for  all  five  of  the  original  Duntley  equa¬ 
tion  coefficients  not  just  the  direct  flux  extinction  coefficient. 


Initially,  asymmetrical  scattering  properties  of  the  individual  leaves 
are  considered.  The  analyses  presented  are  of  theoretical  interest 
since  many  of  the  developments  presented  are  similar  to  the  types  of 
calculations  which  must  be  considered  in  applying  classical  radiative 
transfer  calculations  to  vegetation  terrain  elements,  e.g. ,  to  calculate 
appropriate  canopy  phase  functions.  However,  the  complete  analysis 
presented  is  of  limited  practical  application.  Simpler  approximations, 
e.g..  Chapter  2,  appear  more  relevant. 

The  approach  is  somewhat  complex.  The  general  plan  of  attack  is  as 
follows : 

First,  the  solution  to  the  Duntley  equations  for  arbitrary 
coefficients  is  given. 

Notation  and  definition  of  the  coordinate  system  used  is  then 
given.  The  reader  should  note  that  the  direct  flux  extinction  coeffi¬ 
cient  of  Chapter  2,  k,  is  now  represented  by  and  leaf  slope 
distribution  by  S(0^)  rather  than  f(0^). 

Next,  the  change  in  the  irradiance  in  the  horizontal  plane  caused 
by  scattering  of  flux  in  an  infinitestimal  plane  is  calculated.  This  is 
accomplished  by  dividing  the  coordinate  system  into  four  sections,  so 
that  the  calculation  of  upwelling  and  downwelling  power  emerging  from  a 
single  component  becomes  easier.  Then  an  average  over  all  leaf  inclina¬ 
tions  is  performed  and  the  upgoing  and  downgoing  power  emerging  from  a 
volume  element  calculated  from  this  power  intercepted  by  both  sides  of 
the  component  due  to  downwelling  and  upwelling  flux  is  calculated. 

The  expressions  derived  above  for  the  change  in  the  irradiance  in 
the  horizontal  plane  caused  by  scattering  of  upwelli/i^  and  dowwelling 
flux  are  then  used  to  derive  expressions  for  Then  a 
similar  analysis  used  to  derive  ^23' 
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To  evaluate  “33  >  “se  the  fact  that  the  power  intercepted  by  a 
component  is  proportional  to  the  projection  of  the  component  in  the 
direction  of  the  beam.  So  we  evaluate  the  average  power  intercepted  by 
a  component. 

An  outline  of  a  solution  to  the  problem,  by  using  the  classical 
radiative  transfer  equation  and  an  iterative  procedure  to  find  the 
radiance  emerging  from  the  canopy  at  a  certain  direction  is  then  given. 
However,  this  approach  is  shown  to  be  very  cumbersome. 

Thus,  instead  certain  assumptions  that  simplify  the  model  are  made. 

More  specifically,  that  the  canopy  components  are  assumed  uniformly 

diffuse  Lambertian  scatterers  and  the  upper  reflectance  (transmittance) 

of  the  component  is  assumed  equal  to  the  lower  reflectance  (trasmit- 

tance).  Using  these  assumptions  the  formulas  of  earlier  sections  are 

evaluated  to  find  the  coefficients  a... 

ij 

Finally,  the  bidirectional  reflectance  distribution  of  the  canopy 
is  derived. 


The  Duntley  Equations  and  Their  Solution 

The  Duntley  system  of  differential  equations 

=  -  ajjE(+d,x)  +  aj2E(-d,x)  +  aj3E^(x) 

=  -  a2jE(+d,x)  +  a22E(-d,x)  -  a23E^(x)  (1) 

dE  (x) 

-ir-  =  “33^s^*^ 

Subject  to  the  boundary  conditions: 

E(-d,0)  =  0 

E(+d,l)  =  p  {E(-d,-l)  +  E  (-1)},  where  p  is  the  soil 
s  s  s 


reflectance 
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has  the  solutions: 

Of  X 

E(+d,x)  =  Ae"“  +  Be"®*  +  CE  (0)  e 

s 

Of  X 

E(-d,x)  =  NjAe®*  +  N^Be"®*  +  DEg(O)  e  (2) 

Eg(x)  =  E^(0)  e 

where  a  =  0.5  (022  -  and  a  =  0. 5(0^2  “  “21^ 


_  _  /  -.2  _2  .T  _  Of  +  ni  „  _  01  -  m 

m  -  V  a  -  a  ,  .  N2  =  -5- 


C  = 


“23^^  -  “13  ^“33 
2  „2 
®  -  “33 


D  = 


“13°  ^  “23(“33  ^ 


2  2 
®  "  “33 


and  A,  B  are  found  by 


where  e(0)  =  E(-d,0)  -  D  E^  (0) 


e(s)  =  [p^  (D  +  1)  -  C]  E^  (0)  e"®33 

s 


M  = 


N, 


(1  -  P3  Nj)  e 


-m 


(1  -  P,  N2)  e' 


m 


Establishing  Notation,  Definitions,  and  Coordinate  System 

We  consider  a  single  canopy  element  and  its  interaction  with 
upwelling  and  downwelling  uniformly  diffuse  field.  Let  denote  the 

upwelling  uniformly  diffuse  radiance  of  the  field,  L_  the  downwelling 
uniformly  diffuse  radiance,  E^  (resp.  E_)  the  upwelling  (resp.  down¬ 
welling)  irradiance  evaluated  on  horizontal  plane,  I  the  reduced 

8 
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incident  intensity,  and  E  the  irradiance  resulting  from  reduced 

s 

incident  flow. 


Now  by  definition  n  L±  =  E±  and  E  =  I  cos0  'Ml  ,  where  Ml 
''  s  s  s  s’  s 

is  the  solid  angle  under  which  the  source  of  specular  flux  is  viewed  and 
8^  the  zenith  angle  of  the  source  of  specular  intensity. 

The  direction  K  of  the  incoming  radiance  can  be  expressed  in 
terms  of  the  polar  angles  6  and  4i  as  K  =  (fK|  sin6cos()),  |k| 

A,  ^ 

sin6sin(|),  |  K  {  cos6)  and  if  K  is  considered  to  be  a  unit  vector  (i.e., 
1K\  =  1)  K  =  (sin6cos(t>,sin6sin(t),cos6) .  By  dA  denote  the  vector  down 
the  line  of  maximum  slope  of  the  component  of  area  dA,  and  by  n 


denote  the  unit  vector  perpendicular  to  dA.  I  (K)  indicates  the 
intensity  reflected  from  upper  surface  of  dA,  I  (K)  the  intensity 
reflected  from  lower  surface  of  dA,  I  (K)  the  intensity  transmitted 
outward  from  the  lower  surface  and  incident  on  the  upper  surface,  and 
I  (K)  indicates  the  intensity  transmitted  outward  from  upper  surface 
and  incident  on  the  lower  surface. 

The  bidirectional  transmisstance  and  reflectance  functions  are 


defined  by: 


I  (K):  =  ^  J.  I.„^(K)p.(K,K')  •  K'  •  a  dfi' 

I,  (K):  =  i  J.  I.  (K')  ’  T.(K,K')  •  K'  •  a 
T.  n  ■'i  me  1  ’ 


where  i  =  and  integration  is  performed  on  the  upper  or  lower 
hemisphere  as  defined  by  the  plane  of  dA,  K  is  the  vector  of  incidence 
of  light,  K'  the  vector  of  exit  of  the  light,  and  reflec¬ 
tance  and  trasmittance  functions.  The  trasmittance  function  of  the 
upper  surface  is  related  to  that  of  the  lower  surface  by; 


T.  (K,K')  =  T-  (-K',-K) 
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To  simplify  calculation  we  choose  the  coordinate  frame  as  follows:  the 
z-axis  is  down  the  line  of  maximum  slope  of  the  canopy  component,  the 
y-axis  is  perpendicular  to  dA  (i.e.,  n  =  e^).  Thus,  if  6  is  the 
zenith  angle  of  the  incoming  radiance  vector  and  (|>  its  azimuth  angle, 
then  n  •  K  =  sinBsinij),  since  in  our  coordinate  frame  n  =  (0,1,0). 

It  should  also  be  noted  that  in  this  chapter  leaf  angle  is  denoted 
by  6^  and  leaf  slope  distribution  by  S(6^). 

Calculation  of  Total  Upwelling  and  Downwelling  Power 

Denote  the  inclination  of  the  canopy  element  by  6^.  In  our 
coordinate  system,  partition  space  into  Section  A,  B,  C,  D  as  follows: 

A  =  {(e,(|>)/0  <  6  <  n  -  0^,  n  <  (|)  <  271} 

B  =  {(0,({i)/O  <  0  <  0^,  0  <  $  <  7t} 

C  =  {(0,(t»)/7l  -  0^  <  0  <  7t,  7t  <  <|)  <  271} 

D  =  {(0,<t>)/0j^  <0<7i,  0<(t><7t} 

Then  all  the  upwelling  intensity  emerging  from  the  upper  surface  is 
contained  in  region  A.  All  the  downwelling  intensity  emerging  from  the 
upper  surface  is  contained  in  region  C.  All  upwelling  intensity  emerg¬ 
ing  from  a  lower  surface  is  contained  in  region  B,  and  all  the  downwell¬ 
ing  intensity  emering  from  lower  surface  is  contained  in  region  D. 
Therefore,  all  L_  (resp.  L^)  incident  on  the  upper  surface  is  con¬ 
tained  in  region  A  (resp.  C),  and  all  (resp.  incident  on  the 

lower  surface  is  contained  in  region  B  (resp.  D) . 

Thus,  the  scattered  intensity  caused  by  downward  diffuse  intensity 


L  is  given  by: 
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L 

^  =  r  /a  •  n  •  k'  dn'  (5) 

L 

(e,<>)  =  ^  /a  n  •  K'  dn'  (6) 

L 

ly^g  n  •  K'  dn'  (7) 


L 

It^  (e,«i.)  =  Jg  Xjj(e,0;e ',(!.')  £  •  k'  dn'  (8) 

where  (6,<ti)  is  the  outgoing  direction  of  I,  and  (6 ',<]»')  is  the  incoming 
direction  of  L. 

Let  I  be  the  portion  of  I  ,  I  ,  I  ,  I  directed  upward, 

and  I  be  the  portion  of  I  ,  »  I,.  »  I_  directed  downard.  Then 

^  *£ 

the  power  directed  upward  is: 

P+  =  dA  •  n  •  K  dn  (9) 

and  the  power  P_  directed  downward  is: 

P  =  dA  I  •  a  •  K  dn  (10) 

-  o  o 


The  upwelling  intensity  emerging  from  the  upper  surface  is  equal  to 
I  +  I  and  is  contained  in  region  A.  The  upwelling  intensity  emerg- 
ing  from  the  lower  surface  is  I  '*’1  and  is  contained  in  region  B. 
Similarly  for  downwelling  intensity.  Therefore  (9)  and  (10)  become: 


p+  =  dA  (I  (0,<»  +  (9,4.))  n  • 


K  dn 


+  dA  ;  (I  (0,0)  +  I  (0,0))  a  •  K  dn 

®  *0  X., 


(11) 
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P  =  dA  L  (I  (e,(l>)  +  I  (e,(t>))  n  •  K  dft 

+  dA  L  (I  (0,(1))  +  n  •  K  df2  (12) 

^  '£ 

Let  S(6  )  be  the  distribution  function  of  the  inclination  of  the 
n 

canopy  components,  where  S(6^)  is  normalized  by: 

§  8(0  )  d0  =  1 

n  ■'o  n  n 

Then  the  average  of  upward  and  downward  power  is: 

=  H  '‘»n 

<P.>  '  S  P.(e„)  S(O^)  de„  C14) 

Consider  a  volume  element  dV,  and  an  element  of  horizontal  area 
da.  Then  if  x  is  thickness  of  dV,  dV  =  da  dx.  The  number  of  ele¬ 
ments  in  dV  is  ndV  =  n  da  dx,  where  n  is  the  number  of  canopy  ele¬ 

ments  per  unit  volume.  Thus,  the  upgoing  power  emerging  from  dV  is 
dP^  =  n  da  dx  <P^>  and  the  downgoing  power  emerging  from  a  volume  ele¬ 
ment  dV  is  dP_  =  n  da  dx  <P_>.  Since  the  irradiance  is  defined  as 
power  over  unit  area,  we  have: 

dP  dP_ 

dE.  =  -j —  =  n  <P.>dx  and  dE  =  3 —  =  n  <P  >dx 
+  da  +  -  da 

Hence,  the  change  in  the  irradiance  in  the  horizontal  plane  caused  by 
scattering  of  L_  in  a  layer  of  thickness  dx  is: 

dE  dE_ 

■2 —  =  n  <P  >  and  ^ —  =  n  <P  > 

dx  +  dx 

Analogously  for  the  upward  diffuse  field  we  get: 

L. 

I'  (».♦)  =  JT  -^C  n  •  K  dft 


(15) 
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i;^(e,o)  =  ^  /j)  a  •  k  dn 


i;  (0,<t»)  =  ^  /c  Ty(0,O;e',4.')  a  •  K  dft 
and  using  (15)  -  (18)  we  get  that  the  power  is: 

p;  =  dA  /.  (I'  (0,(|.)  +  i;  (0,<|.)i  a  •  K  dn 

*  'o 

+  dA  /  ir  (0,«)  +  i;  (0,4.)]  a  •  K  dn 

o 

p_'  =  dA  j  ir  (0,4.)  +  i:  (0,4.)]  a  •  K  dn 

>u 

+  dA  /  (I'  (0,4»)  +  i:  (0,4»)]  a  •  K  dn 
“  *£ 

and  averaging  over  0^  we  get: 

<P>  =  i  K  (®  )  S(e  )  d0 

+  n-'o  +  n  n‘'  n 


<?'>  =  §  p'  (0  )  s(0  )  d0 

n-'o  -  n  n  n 


Thus, 


dE^  dE_' 

—  =  n  <P'>  and  -j —  =  n  <P'> 
dx  +  dx  - 


where  the  primes  indicate  irradiance  caused  by  L^.  We  now  calculate 
the  total  power  incident  on  dA.  By  definition  the  incident  power  P. 


is  given  by: 


Pine  =  fi  (®-<'5  n  •  K  dn 
inc 


The  downwelling  diffuse  radiance  L_  is  nonzero  only  in  the  region  A 
and  B,  and  since  it  is  uniformly  diffuse  it  has  to  be  a  constant.  Also 
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dfi  -  sin0d6d(|>,  therefore  the  total  power  delivered  to  the  upper  surface 
of  dA  by  L_  is: 

Py  =  dA  •  L_  n  •  K  dft  =  dA  ♦  L_  sin  6sin())d6d(|) 

since  the  dowelling  diffuse  radiance  L_,  incident  on  the  upper  side  of 
dA,  is  contained  in  region  A. 

Evaluating  the  integral  for  we  get: 

n-e 

Py  =  dA  •  L_  sin  6sin(|)d6d(|»  =  (n  -  0^  +  %sin(26^)]  dA  •  L_  (23) 

In  the  same  way  the  total  power  P^  delivered  by  L_  onto  the  lower 
surface  is: 

P^  =  dA  •  L_  /g  sin^6sin4)d6d(|i  =  (6^  -  ^sin(26^)]  dA  •  L_  (24) 

the  total  power  Pj  delivered  by  onto  the  upper  surface  is: 

P'  =  dA  •  L  sin^6sin()id0d(|>  =  [0  -  ^sin(20  )]  dA  ♦  L.  (25) 

and  the  total  power  Pj  delivered  by  1^  one  the  lower  surface  is: 

Pj^  =  dA  •  /g  sin^0sin<|id0d(^  “  ”  ®n  ^  %sin(20^)]  dA  •  L  (26) 


Hence,  the  total  power  intercepted  by  both  sides  from  L_  is  ~ 
PjU  =  7lL_dA  and  the  total  power  from  is: 


=  nL_j_dA 


Calculation  of 

The  portion  of  lost  over  dx  because  of  either  absorption  or 
redirection  downward  of  is  equal  to  Therefore  i® 
the  difference  of  the  total  upwelling  power  intercepted  by  dA  minus 
the  amount  of  power  emerging  upward  from  dA. 
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Similarly,  ^22^-  portion  of  £_  lost  over  thickness  dx 


due 

to 

absorption  or  redirection 

of 

E.. 

“12®- 

is  the  portion  added  to 

over  dx 

due 

to  scattering 

of 

E_ 

into 

E^,  and  <*21®+ 

change 

in  E_ 

due 

to  scattering 

of 

into 

E_.  So  using  the  pre- 

viously  derived  formulas  for  total  power  we  get: 
=  n(L^ndA  -  <P^>) 


^22!^-  ~  -  <P_>) 


=  n<P_'> 
dE^ 

"12®-  =  dT  =  "  "V 


and  since  E^.  =  nL^  we  get: 


11 

dA 

n 

<P^> 

12 

“  1 

n 

n 

'21 

K  " 

'22 

=  n 

TdA 

<P 

L 


+  " 


(27) 

(28) 

(29) 

(30) 

where  <P+>,  ate  given  in  formulas  (13),  (14),  (21),  and  (22). 

The  Derivation  of  “23 

The  change  in  upwelling  diffuse  intensity  over  dx  due  to 
scattering  of  specular  intensity  is  aj2E(8,o),  where  E(s,o)  is  the 
specular  flux  at  the  top  of  the  canopy,  and  is  given  by: 


E(s,o)  =  I(s,o)  cos© 


/in 


In  the  same  way,  a22E(s,o)  is  the  change  in  downwelling  diffuse  inten¬ 
sity  over  dx  due  to  scattering  of  specular  intensity. 
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The  scattered  intensity  from  volume  element  AV  on  which  intensity 
I(s,o)  is  incident  is  given  by: 

I  (cos0,(Ji)  =  \n  P(cose,<t>;  cos6  ,(J»  )  I(s,o)  ~  Afi 

S  S  S  ijA  S 

where  P(cos6,<|>;  cos6^,(^^)  is  the  phase  function. 

The  power  scattered  from  volume  element  AV  is: 


P^  =  AA  I^(cos0,<ji)  dQ 

=  k,n  AV  I(s,o)  Afl  P(cos0,((»;  cos0  ,({»  )  dfi 
S  s  s 


where  +  means  upwelling  and  -  downwelling  power  and  +h  denoted 
integration  over  the  uppwer  hemisphere,  -h  over  the  lower  hemisphere. 
Now  since: 


AV 


and 


E(s,o) 

COS0 

s 


I(s,o)  •  Afl^ 


we  have: 


Thus: 


^  =  E(s.o)  r 

dx  A71COS0  -^ih 
s 


P(cos0,(t»; 


COS0  ,(|)  )dfl  =  a.»E(s,o) 

S  S  X 


“l3  ~  Anp^  •^+h 
“23  “  A^  -^-h 

wh^te  M  =  COS0,  =  cosO^. 

Evaluation  of  0^3 

The  power  intercepted  by  a  component  is  proportional  to  the 
projection  of  the  component  in  the  direction  of  the  beam. 

Let  X  be  a  vector  of  magnitude  equal  to  the  area  dA  of  a  compo¬ 
nent  and  direction  perpendicular  to  the  plane  of  the  component.  Let  ^ 
be  the  unit  vector  in  the  direction  of  the  incoming  beam. 


In  a  cartesian  x,y,z  system  with  z  vertical  we  have: 

^  =  Z  A.  £.,  i=x,y,z,  ^  •  A  =  A^  =  dA 
i 

^  =  Z  R.  e.,  i=x,y,z,  ^  S  =  1 
i 

^  ^  ^ 

With  8  ,  &  ,  e  unit  vectors  on  x,  y,  and  z  axis  respectively. 

X  y  z 

In  polar  coordinates  we  have: 

X  =  tA\(sin6  cosA  e  +  sin0  sin*  £  +  cos6  £  ) 

'  »  n  ^n  X  n  ^n  y  n  z 

^  =  sin6  cosd)  £  +  sinS  sin*  £  +  cos0  £ 

s  ^s  X  s  ^s  y  s  z 

Thus: 


I  ^  =  I  A|  [  jsin0^sin0^cos((|>^-(^^)  +  cos0^cos0^|] 

By  basic  algebra  we  see  that  if  |^  *  ^|=  |A(  {|F((t)^-(|)^)| }  with 

F((li  -d>  )  =  sin0  sin0  cos (4  -4  )  +  cos0  cos0  then  for  4  ”4  ^  6  : 
^n  ^s  ns  ^n  ^s"^  ns  ^n  ^s  o 

F(4  "4  )  =  0  and  for  4  “4  >5  :F(4  -4)^0  where: 

rn,  for  0  <  0  +0  S  5 
\  ’  n  s  2 


rn,  for  0  ^ 

°  1  cos  ^  (-cotO^ 


cot0  ] ,  for  0=0  +6  <  n 

s’  2  n  s 


Hence,  if  we  assume  that  the  distribution  functio  for  4  is  1^(4  )  =  1 

n  n 

(i.e.,  assume  azimuthal  symmetry),  we  can  calculate  the  average  of 
X  •  ^  as  follows: 


<  ^  '  S  >  =  -  {/  °  [sinB  sin0  cos (4  “4  )  ■•■ 

4  n  o  n  s  n  s 


+  cos0^cos0g]  <l(4jj"4g)  I-sin0^sin0^cos(4„”4e)  “ 


n  s  'n  's 


COS0  CO80  ]  d(4  “4  )}  = 
n  s  ^n  ^s 


=  {sinO  sin0  8in6  +  (6  -  ?)  cosO  cos0  }  S(9  )  d0 

n‘  n  s  o  o2  n  s’  n  n 
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Averaging  next  over  6^,  we  find  that  the  average  power  for  component  is: 
.  .  An  "/2-0 

<  A  •  R  >-  .  =  — r  |A|  cos6  /  ®  S(0  )cos0  d0  + 

0  2  '  2  s  •'o  n  n  n 

n’^n  71 

+  sin0  a  sin0  sin6  S(0  )d0  + 

s  •'n/2-0  n  o  n  n 

s 

+  COS0  a  “  ?)  COS0  S(0  )d0  } 

s  ■*n/2-Q  '  o  2  n  '  n  n' 

s 

where  6  =  cos  ^[-cot0  cot0  ].  Finally, 

o  ns  ^  * 


33 


n  •  •  R>  0  ,(j) 

_ n*^n 

COS0 


(33) 


Evaluation  of  Emergent  Intensity  in  the  General  Case  of  Plane  Parallel 
Problems 


Let  p  =  COS0,  Then  the  radiation  transfer  equation  becomes: 
M  -  F(t,p,(|)) 

where  the  optical  thickness: 

T  =  kpdt 


k  =  scattering  coefficient  and  p  =  density  of  material.  If  the  optical 
thickness  is  finite  then  the  solution  becomes: 


I(X,M,<|))  =  I(Tj,M,(|»)e 


-(Xj-x)/M  Xj 


+  F(t,p,(|.)e"^'^'^^/^'  ^ 


In  particular  for  emergent  intensity  we  have: 

.t/u  dt 

L(m,<|»)  =  I(o,M,<|>)  =  +  /o  ®  F(t,M,<|>)  ^  (34) 

where  the  source  function  F  is: 

F(0,4»)  =  hjl  f]  (35) 


To  solve  (34)  we  use  the  following  iterative  approach: 


First  we  make  an  initial  guess  for  I,  by  using 
the  Duntley  equations.  Then  given  this  initial 
guess  for  I,  we  can  calculate  F  from  (35)  and 
then  update  I . 

Now  let  fhe  view  angle,  the  direction  of  coming  specular 

flux,  =  the  soil  reflectance.  Assume  components  in  canopy  are 

Labmertian  scatters.  Then  the  upward  radiance  L  of  the  canopy  for  the 
view  angle  (|J,<|i)  will  be: 

1  -t/u 

L(|J,(j>)  =  I(o,M,(t))  =  I(Xj,M.<l))e  J  -To  F(t,p,<}.)e  dt 

Now  I(T^,p,<ti)  is  the  radiance  from  the  soil  that  goes  through  canopy 
and 

I(Tj,|J,<|l)  =  Yg  I(+d,Tj) 

To  find  the  initial  I(o,|j,(ti)  we  assume  that  the  medium  (canopy)  is  a 
conservative  case  of  perfect  scattering,  i.e.,  all  incident  radiation 
becomes  diffuse  flux;  there  is  no  absorption.  In  that  case  the  upward 
radiance  will  be: 

I(T,Mj,(t.j)  =  I(+d,x)  +  I(-d,x)  +  I(s,|J,,,<l>^) 

But: 

E(±d,x)  =  nL(±d,x)  =  n  I(±d,x) 

so: 

I(+d,x)  =  I(-d,x)  = 


and  also: 


I(s,x,|J^,(|»^)  =  I(s,o,M^,<|»Q)e 


■(T-Tj)/M„ 


E^(o)  -(t-Xj)/M„ 

-  e 


Thus : 
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rl  r2n 


}  d(|»'  djj 


E  E 

=  ^  /q"  P(|J,(tl,^J ',<!>')  +  sli  Jq"  P(M,0,M'.4'')d<Ii'd|j' 

An  An 


,1  ,2n 

Anjj  ® 


X.l  xr  P(M.<l>.M%<l>')d4>'  d)J- 


And  so: 

F(e,(i>)  =  F(t,M,<t>)  =  kn  /I  ;2"  P(^J,,|,,^J',,|,')d^J'  d<i.' 

+  E(+^  jO^  j2n  d(|.'} 

-  e  P(M,<t),M^,V 

Now  if  we  know  the  phase  function  we  can  find  I(o,p,<t>)  =  L(p,(l>). 


Development  of  a  Simplified  Model  for  Lambertian  Canopy  Components 

To  simplify  the  model  we  assume  first  that  the  canopy  components 

are  uniformly  diffuse  Lambertian  scatterers  (i.e,  that  the  bidirectional 

scattering  functions  are  constants). 

Let  G  =  0  -  ^sin(20  ).  Then: 

n  n 

sin0sin(|idQ  =  n  -  0^  +  %sin(20^)  =  n  -  G 

/_  sin0sin(|idQ  =  0  -  %sin(20  )  =  G 

D  on 

sin0sin(|>df}  “  ”  ^sin(20^)  =  G 

Jjj  sin0sin(|idQ  =  n  -  G 

Now  since  n  •  K  =  sin6sin(^,  using  relations  (5)  -  (8)  for  L_,  we  get: 

“  JT  Xa  Pp  "  •  Kdn  =  sin  0sin<|»d0d(|»  =  (n  -  G) 
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L  T 

I  (6,(1))  =  (n  -  G) 

'"u 

^_Po 

I  iQA)  =  ^  G 
*£ 

« 


Similarly  using  relation  (15)  -  (18)  for  L^,  we  get: 

-  ^+Po 

1 '  =  —  p  n  •  K  dn  =  — ^  G 

y  71  •'C  *^0  71 

'u 


Lx 

i:  (6,(1))  =  ^  G 

*■0 


i:  (QA)  =  ^  (rt  -  G) 

*£ 


K'^SL 

i;  (6,(1))  =  (n  -  G) 


Using  the  above  equations  in  relations  (11)  and  (12)  we  get  that  the 
power  due  to  downwelling  intensity  is: 


dAL 


[p^j(7r  -  G)^  +  x^  G(71  -  G)  +  Pj^G  +  in  -  G)G] 


dAL  2  o 

P_  =  [Py(7i  -  G)G  +  X^G  +  Pj^G(7T  -  G)  +  x^j  (n  -  G)  ] 

Similarly  the  power  due  to  upwelling  intensity  is: 

dAL ,  y  y 

p;  =  [p^jG(7T  -  G)  +  X^(7l  -  G)^  +  P^(7t  -  G)G  +  X^G  1 

dAL .  y  y 

P;  =  (p^G^  +  ij^(n  -  G)G  +  Pj^(n  -  G)  +  T^(7i  -  G)G] 
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where  ®  ~  ^sin(26^).  Thus  averaging  over  6^  and  using  relations 

(13)  and  (14),  we  get: 

<P+>  =  [D(p^  -  +  Pa  -  t^)  +  R(-2p^  +  +  t^)  +  pj 

<P->  =  [D(-p^-Pjg  +  +  l^)  +  R(p^  +  p^  -  2T^)  +  Ty] 

<P;>  =  TidAL^  [D(-P^-Pa  +  +  t^)  +  R(p^  +  p^  -  2x^)  +  x^] 

<P'>  =  ndAL  [D(p  +p  -I  -t)+  R(t  +  t  -  2p  )  +  p  ] 

+  021LU  u££2 


where: 

R  =  \  t®n  “  ^sin(26^)]  S(0^)  de^ 

71 

D  =  ^  [6  -  ^sin(0  )]^  S(0  )  d0 

^3  -'o  ‘  n  n  n  n 

Now,  relations  (27)  -  (30)  give: 

<p;> 

“ll  =  “  t<lA  -  =  ndA  [1  -  (Vt^-p^-p^)D  -  (p^+p^-2Tj^)R  -  t^] 


a 


12 


a 


21 


n<P^> 

=  ndA  l(p^  +  “i  ♦  +  (!„  *  t. 


n<P  > 

sir  = *  pj  ■  -  ’u>“  *  <’u  *  h 


■  2Pu>*  *  fu> 

-  2Pj^)R  +  Pj|,l 


<P  > 

«22  =  “  IdA  -  =  ndA  [1  -  Xy  -  (x^+x^-  Pjj-p^)D  -  (P/Pa-2x^j)R1 

Assuming  that  p^  =  and  x^  =  Xj^  we  get: 

fljj  =  ndA  [1  +  2(x-p)D  +  2(p-x)R  -  x]  =  022 

0^2  =  ndA  [2(p-x)D  +  2(T-p)R  +  p)  =  (*21 


/ 


with 


^  G(0„)  S(0^)  d0^ 

_z  o  n  n  n 

71 

D  =  ^  G^(6„)  S(e„)  d0„ 

_3  o  n  n  a 

71 


G(0  )  =  0  -  ^siii(20  ) 

n  n  n 

Evaluation  of  023  for  case  =  p^,  =  Tj^ 

The  intensity  due  to  reflection  from  upper  surface  is: 


I  ^  p.  .s 
*u 


(37) 


The  intensity  due  to  transmission  from  the  upper  surface  to  the  lower 
surface  is: 


I  An  T  . 
o  o  u 


fK  •  n\ 


(38) 


The  intensity  due  to  reflection  from  lower  surface  is: 


I  An  p. 
o  o 


U  *  n( 


(39) 


The  intensity  due  to  transmission  from  lower  to  upper  surface  is: 


I 


I  An 
o  o 


71 


lie  ■  ai 


(40) 


where  I  K  •  n  I  =  cosy  =  sin0  sin0  cos(di  -(!>)+  cos0  cos0  and  v 
<  I  *0  s  n  ^s  ^n  s  n  *0 

A 

is  the  angle  between  and  h. 

The  conditions  for  upper  and  lower  incidence  can  be  seen  to  be  the 


following: 
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cosY^  is  positive  if  ^  ^  ^  ~  ®s  ®  ^  *^n  ^ 

or  if  (5-6  <6  <5  and  (|>  -  6  <()>  <4i  +  6) 

2  s  n  2  o  ^n  o 

Also  COSY  <0  if  (?“0  <0  <?  and  ())  +  6  <(|»  <4  +2n-6) 

*o  '2  s  n  2  ^s  o  ^n  ^s  o 


where 


6  =  cos  (-cot0  cot0  ) 

o  s  n 


The  upwelling  power  reflected  from  upper  side  is 


P  =  I  dA  /.  sin0  sin<t>  d£l 
*u  *u  * 


Similarly 


P  =  I  dA  /-  sin0  sin(|)  dQ 
»u  ^ 


P  =  I  dA  sin0  sin<j>  dfl 


(41) 


P_  =  I_  dA  /_  sin0  sin<l>  dQ 


Using  (37)  and  (Al)  we  get: 


I  Aft  cosy  dA 
o  o  *o 


P.  (0  ,<f  )  = 

u  n’^n 


[p  n  +  (t  -  p  )  G] 
u  u 


E  (0)  cosy 

-£-j_  __  dA  Ip^A  ♦  (,„  -  p^)  Cl 


I  An  cosy  dA 

-  'V"  *  %  ■ 


E  (0)  cosy 

“  I'u"  *  ‘Pu  -  'u>  '=> 


Similarly  for  elements  illuminated  on  the  lower  side  we  have: 
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.  I  ^  cosy  dA 

^  IV  "  (P£  ■ 


E^  (0)  cosy^  dA 
n  cos6 


I V  ^ 


I  A^2  cosy  dA 

h  °  °  ^  -  IPo"  +  -  Po)  G] 


71 


E^  (0)  cosy^  dA 

71  COS0 


IP£”  ^^£  ■  P£^ 


Thus  the  total  upwelllng  power  by  all  elements  in  AV  is  given  by: 

p+  =  S|V  j  ^  s  j.271  j  j  P+  ^  j  j  Q 

2  •’o  ■'o  n  '^n  u  n’^n  ^n  n 

7t 

"  <l'/U  4’-^  <  («„-♦.>  "♦n  ■“’„ 

s  o 

/»  <|i  +27t-6  . 

V*6  °  ^‘VXV  fi!  «„■♦„>  "♦n  <«) 


s  s  o 


and  the  total  downwelling  power  is  given  by: 

...  71/2-0  - 

p-  _  n^  J.  S  J^Tt  J  ^  J  ^  ^ 

„2  •'o  •'o  '  n"^  '^n'  u  n’^n"^  ^n  n 

A 

<l>.+6 


*  V-5°  “‘♦n^  ‘’u  '*'»«  •‘®n 

s  s  o 
(t>  +271-6 

"  •^/2-e  -^Z+a  °  ‘'♦n 


(43) 


Now  as  before  023  ’  ^ 

s  s 


Assuming  that  the  distribution  function  D(0^)of  the  azimuths  of  the 
leaves  is  equal  to  1  we  get: 
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=  —qr -  ^  S(0  )(np  +  (t  -  p)  G]  (27t  cos0  cos0  )  d0 

COS0  on  s  n  n 

s 

+  jJJy2-0  I^P'*'(T“P)GJ  (2sin0^sin0^sin6^+26^ 

+  a  S(0  )[7tT+(p-T)G](2(n-6  )  COS0  COS0  -2sin0  sin0  sin6  )d0 

7t//-o  n  osnsno  n 


Also 


a 


ndA 


71/2-0 


23  „3  -  'o 

71  cos6 

s 


S(0  )[nx  +  (p  -  t)  G]  271  COS0  COS0  d0 
n  vr-  ^  j  g  n  n 


+  /!!/o  a  S(0  )[7iT+(p-T)G](2sin0  sin0  sin6  +  26  cos0  cos0  ) 
Ti/z-tj  n  s  n  o  o  s  n 


d0 


+  /|J/2-0  ”  2sin0^ 


Evaluation  of 

In  the  same  way  as  in  previous  section  we  find  that  = 

n*<  >Q  X  /  Q  where 
0  ,<b  /cos0 
n  ’  ^n  s 


"0  ,4.  =  i  "°®®S  J”. 

n’  n  ^ 


71/2-0 

®  S(0  )  COS0  d0 
o  n  n  n 


fji/i 


+  sin0  S-J/n  a  sin0  sin6  S(0  )  d0  +  cos6 
s  •'71/2-0  n  o  n  n  s 

s 

-  i)  '»=»n 


and  6  =  cos  (-cot0  cot0  ). 

o  ns 


Calculation  of  the  Bidirectional  Reflectance  Distribution  Function  of 
the  Canopy 

Due  to  the  difficulties  involved  in  using  the  phase  functions  to 
calculate  the  bidirectional  reflectance,  we  will  use  here  an  approach 
similar  to  that  used  in  the  SRVC  Model  (Smith  and  Oliver,  1972). 


f 


-63- 


We  have  that  the  reflected  and  transmitted  radiation  in  direction 
6  due  to  incident  flux  is:  (see  Goudriaan,  1973,  p.  A8). 

(i)  Due  to  specular  flux  E  ,  the  radiance  is: 

s 

4 

Lp  =  E„  {(p  +  t)  — ^  [a,  (2a  -  na^  -  na')  +  a,  sina'  (2a„  -  n) 

tg  a  ICC  c  c  z  cc 

+  a_  sina  (2a'  -  n)  +  2a,  sina  sina']  +  2a- p}  sin6  (44) 

J  c  c  4  c  c  1 

(ii)  Due  to  downwelling  irradiance  E(-d): 

h(-d)  =  ^  A  S  f“l  ^2“c  “'c  -  "“c  -  "“c^ 

271 

+  a^  sina^  (2a^  -  n)  +  a^  sina^  (2a^  -  n)  +  2a^  sina^  sina']  dS 
+  2p  I  “l  <>61  Sine  (45) 

(iii)  Due  to  upwelling  irradiance  E(+d): 

It<M)  =  {(p  ♦  t)  2  |»,  (2«^ 

2K 

+  sina^  (2a^  -  7t)  + 

+  2t  ^  tA6) 

2 

where,  a,  =  cos6  cosO  cos  6 

'  I  O  M 


=  cose 

cose 

cos^e 

S 

n 

=  cose 

sine 

sine 

cose 

s 

n 

) 

=  sine 

cose 

sine 

cose 

s 

n 

1 

=  sine 

sine 

•  2„ 
sin  e 

cos  ^  (-cote  cote  ),  if  ?  -  e  <  e 
n  s’  n  n  8 


71,  if  63 


(47) 
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a' 

c 


f  cos  ^  (-cot0  cot6),  if  6  >  -  -  0 
\  n  ’  71  n 

/  Tt.  if  e  s  5  -  e^. 


Thus  the  total  radiance  emitted  from  a  component  of  inclination  0^  will 
be: 

“  =  %  *  ^(td)  ♦  ■■E(-d)  *  r  *  r  *  r 

where  u  ,  v, ,  w  are  functions  of  0,  0  ,  a.,  a  ,  a'.  For  exact 

111  ’  s  i’  c'  c 

expression  of  u^,  v^,  Wj  see  Appendum  A. 

Now  if  0^  ^  0,  then  71/2  -  0^  ^  rt/2  -  0  and  in  averaging  AL(0^) 

over  0^,  we  need  to  subdivide  the  interval  (0,  n/2)  into  the  subinter¬ 
vals  (0,  n/2  -  0  ),  (n/2  -  0^,  n/2  -  0),  (n/2  -  0,  n/2).  Then  the 

s  s 

values  of  a^,  a'  in  the  expressions  of  u^,  v^,  Wj  become: 

in  w,  (a  ,  a')  we  have: 

1  c  c 

(i)  a  =  n,  a'  =  n,  if  0  §  0  ^  ^  ~  d 

c  ’  c  ’  n  2  s 

(ii)  a  =  cos  ^  (-cot0  cot0  ),a'=n,  if  ?-0  S0  ^?-0 

c  s  n’c  2sn2 

(iii)  a  =  cos  ^  (-cot0  cot0  ),  a'  =  cos  ^  (-cot0  cot0  ), 

c  s  n  ’  c  n  ’ 

if  f  -  0  <  t  <  f 
2  n  2 

in  u,  (a')  we  have: 

1  c 

(i)  a'  =  n,  if  0  ^  0^  g  5  .  0 

(ii)  a'  =  cos  ^  (-cot0  cot0  ),  if  ^-0S0 

c  n  2  n  2 
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in  V 

1  («c) 

we 

have: 

(i) 

0 

a  = 

if  0^0  sH-e 

c 

n  2 

(ii) 

0 

a  = 
c 

cos 

^  (-cot0  cot0^) ,  if  j  - 

‘  0 

Similarly  for 

0  <  0. 
s 

Let 

S(6n) 

be 

the  distribution  function 

of 

the  leaf 

aging 

over 

0 

n 

we  get  that: 

AL(layer)  =  |  dA  w  +  u  +  vl  ,  if  8  ^  6 

n  In  n  n  \  ’  s 

and 

AL(layer)  =  |  w'  +  u  +  v^  ,  if  6^  <  6 

and  the  exact  expressions  of  w,  w',  u,  v  are  given  in  the  Appendum  B. 

Kx 

The  probability  of  achieving  line  of  sight  through  the  canopy  is  e  , 
where 

K  =  I  S(0  )|  X'r  Id0 
n  •'o  n  I  •  n 


where  X  is  a  perpendicular  to  the  leaf  vector  and  r  is  a  unit  vector 
in  the  direction  of  view.  Thus  in  the  same  way  as  for  we  get  that 

n*<  X‘r  >, 


K  = 


0  ,(t) 
n’^n 


COS0 


where 

<  X*r  >Q  .  =  i  dA  [5  COS0  ®  S(0  )  cos0  d0 

0  ,<>  71  2  o  •'o  n  n  n 

n  n 

+  sin0  a  Sind  sin6  S(0  )  d0  +  cos0 

0  ■’71/2-0  n  o  n^  n  o 


Jj/te  <»o  -  5)  '»■»«  '’».l 
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and 

6  =  cos  ^  (-cot0  cot6). 

o  n 

Hence  ALp^  (x)  =  ALp^  (layer)  e^ 

and  so  ALp^  (outside)  =  ALp^  (layer)  e  dx. 


The  solution  to  the  system  of  differential  equations  (1)  is: 

r  tA.A  \  A  a.  t.  “®x  /A  kx 

E  (+d,x)  =Ae  +Be  +Ce 

E  (-d,x)  =  Mj  A  e"“  ♦  B  e"“  »  D  s'" 

E,  (X)  =  E.  (0)  x'“  X 
s  s 

where  A,  B,  C,  D,  Nj ,  and  N2  are  defined  earlier  in  terms  of  j • 
Thus  the  outgoing  irradiance  from  depth  x  is 

AL(layer)  (x)  =  ^  u  +  ~  v  +  w  e*^  (49) 

Kx 

where  u,  v,  and  w  are  given  previously  and  e  is  the  probability 
that  the  light  will  go  through  the  canopy.  To  get  the  total  AL  we 
integrate  with  respect  to  x 

AL(layer)  =  AL(layer)  (x)  =  ^  E  (+d,x)  e^  dx 


+  - A 

n  •'-1 


E  (-d,x)  e*^  <1*  ^  E  (x)  e*^*  dx 

•I  *  1  s 


_  u/a  (1  -  e"^*^  "^  '"^)  B  (1  -  e"^*^  ‘  ”^)  C  (1  -  e"^^  *  ‘'^)] 
n\  K  +  m  K-m  K  +  k  / 


V  ini  ^  ~ 

7T  \  K  + 


e"(K  *  '"^)  N2  B  (1  -  e"^*^  " 


m 


K-m 


)  D  (1  -  e-«  * 

K  +  k  / 


+ 


w 

n 


1 


e-(K  +  k) 

K  +  k 


(50) 


-67- 


Now  the  radiance  contribution  from  the  soil  can  be  calculated  using 
the  hemispherical  reflectance  of  the  soil  and  the  flux  values  of 

E  (-d,-l),  Eg  (-1): 

TtLsoii  =  Pg  (Eg  (-1)  +  E  (-d,-l))  (51) 

Thus  the  total  directional  radiance  of  the  canopy  is 
nL  (0,0g)  =  nL  (layer)  + 

Let  Ej  (0^)  =  E  (-d,o)  +  E^  (o)  cos0^ 

Then  the  directional  reflectance  is: 

R  (0,6^)  =  nL  (0,0g)/Ej  (0^).  (52) 


/ 
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APPENDUM  A 


EXPRESSIONS  FOR  u^,  v^,  and 


“  =  ‘ej  »  *  hi-i, '  r  *  r  *  r 's  ■ 

We  give  below  the  exact  expressions  of  u^,  v^,  and  Wj. 

2  2  2  > 

Uj  =  ( (p  +  ®n  ^  COS0  cos0^  +  sinaj  sin0  sin0^] 

2  >■ 

+  ((2a'  -  7t  COS0  cos  0  +  2sina'  sin0  sin0  cos0  ). 

c  n  c  n  n 

n 

2  “1 

J_  COS0  cos  (-cot0  cot0  )  d0  +  ((2o'  7i)  COS0  COS0  sin0 

n  Q  n  c  n  n 

2  ■  ®n 


n 

2  2 

+  2sina'  sin0  sin  0  )  f 

c  n  "^Tt  . 

2  ■  ®n 


1  .  £0§!6 


sin^0 


n 
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APPENDUM  B 

EXPRESSIONS  FOR  w,  w' ,  u,  v 
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Simulation  of  solar  radiation  absorption  in 
vegetation  canopies 

D.  S.  Kimes  and  J.  A.  Smith 


A  solar  radiation  canopy  absorption  model,  including  multiple  scattering  effects,  was  developed  and  tested 
for  a  lodgepole  pine  (Hnus  contorto)  canopy.  Reflectance  above  the  canopy,  spectral  transmittance  to  the 
ground  layer,  and  geometric  and  spectral  measurements  of  canopy  elements  were  made.  Relatively  large 
differentials  occurred  in  spectral  absorption  by  canopy  layers,  especially  in  the  photosynthetically  active 
region,  as  a  function  of  solar  zenith  angle.  In  addition,  the  proportion  of  total  global  irradiance  absorbed 
by  individual  layers  varied  greatly  as  a  function  of  solar  zenith  angle.  However,  absorption  by  the  entire 
canopy  system  remained  relatively  constant 


I. .  Introduction 

The  manner  in  which  a  vegetation  canopy  absorbs 
solar  radiation  has  an  important  effect  oh  both  the 
thermal  properties  and  the  photosynthetic  efficiency 
of  the  canopy.  Thus,  an  understanding  of  these  prin¬ 
ciples  is  important  in  remote  sensing  with  respect  to 
military,  agricultural,  forestry,  and  ecological  applica¬ 
tions.  In  recent  years  the  thermal  region  of  the  elec¬ 
tromagnetic  spectrum  has  received  keen  interest  in  the 
field  of  remote  sensing.  This  region  may  add  valuable 
additional  information  to  make  inferences  concerning 
the  characteristics  of  vegetation  canopies.  However, 
before  the  thermal  emission  characteristics  of  a  canopy 
can  be  understood,  the  manner  in  which  the  canopy 
'isorbs  solar  radiation  must  be  studied.  In  the  field 
of  agriculture,  there  is  strong  evidence  that  the  solar 
radiation  distribution  within  a  canopy  as  a  function  of 
canopy  structure  strongly  affects  the  productivity  of  the 
canopy.*-® 

Physically  based  mathematical  models  serve  as 
convenient  tools  for  studying  the  complex  radiation- 
vegetation  interactions.  The  objective  of  this  investi¬ 
gation  was  to  develop  a  model  to  study  spectral  and  total 
solar  radiation  absorption  in  vegetation  canopies  as  a 
function  of  solar  zeni^  angle.  The  following  describes 
the  absorption  model  and  the  application  of  the  model 
to  a  lodgepole  pine  (Pinus  contorta)  canopy  at  Lead- 
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ville,  Colorado,  for  which  a  data  base  was  collected 
during  1977.^  A  complete  study  site  description  is  given 
by  Ranson  et  al.*  The  specific  canopy  modeled  con¬ 
sisted  of  a  cluster  of  four  Mgepole  pine  trees  with  these 
mean  statistics:  6.0-m  height,  30-yr  age,  13.2-cm  diam 
breast  height,  and  a  surrounding  stand  of  102-m®/ha 
basal  area. 

II.  Solar  Radiation  Canopy  Models 

Several  deterministic  models  have  been  developed 
to  study  the  interactions  of  solar  radiation  within  veg¬ 
etation  canopies.  Allen  and  Richardson,®  Alderfer  and 
Gates,®  and  Suits’*  have  adapted  a  system  of  simulta¬ 
neous  differential  equations,  develop^  by  Kubelka  and 
Munk,®  in  various  ways  to  vegetation  canopies.  The 
Suits  model  includes  geometric  effects  and  predicts 
non-Lambertian  characteristics  of  vegetation  canopies. 
Chance  and  LeMastei®  derived  a  light  absorption  model 
for  vegetative  plant  canopies  from  the  Suits  reflectance 
model. 

Another  approach  developed  by  Oliver  and  Smith*® 
is  the  solar  radiation  vegetation  canopy  (SRVC)  model. 
This  m  jdel  simulates  the  solar  radiation  flow  through 
the  canopy  by  using  Monte  Carlo  techniques.  This 
stochastic  mc^el  predicted  the  diurnal  apparent  di¬ 
rectional  spectral  reflectance  of  a  vegetation  canopy. 

The  ray  tracing  technique  used  in  the  SRVC  model 
is  advantageous  when  applied  to  solar  radiation  inter¬ 
actions  within  vegetation  canopies  for  several  rea¬ 
sons. 

(1)  The  model  can  accept  any  relevant  parameter, 
such  as  the  number  of  components,  the  optical  prop¬ 
erties  of  the  components,  and  the  spatial  dispersion  and 
density  of  the  components. 

(2)  Diffuse  skylight  is  treated  as  a  set  of  independent 
source  vectors. 
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Thermal  IR  exitance  model  of  a  plant  canopy 

D.  S.  Kimes,  J.  A.  Smith,  and  L.  E.  Link 


A  thermal  IR  exitance  model  of  a  plant  canopy  based  on  a  mathematical  abstraction  of  three  horizontal 
layers  of  vegetation  was  developed.  Canopy  geometry  within  each  layer  is  quantitatively  described  by  the 
foliage  and  branch  orientation  distributions  and  number  density.  Given  this  geometric  information  for  each 
layer  and  the  driving  meteorological  variables,  a  system  of  energy  budget  equations  was  determined  and 
solved  for  average  layer  temperatures.  These  estimated  layer  temperatures,  together  with  the  angular  dis¬ 
tributions  of  radiating  elements,  were  used  to  calculate  the  emitted  thermal  IR  radiation  as  a  function  of 
view  angle  above  the  canopy.  The  model  was  applied  to  a  lodgepole  pine  (Pinu*  eontorla )  canopy  over  a 
diurnal  cycle.  Simulated  vs  measured  radiometric  average  temperatures  of  the  midcanopy  layer  corre¬ 
sponded  within  2*C.  Simulation  results  suggested  that  canopy  geometry  can  significantly  influence  the  ef¬ 
fective  radiant  temperature  recorded  at  varying  sensor  view  angles. 


I.  Introduction 

The  thermal  IR  region  (3-20  ftm)  of  the  electromag¬ 
netic  spectrum  may  provide  valuable  information  about 
the  characteristics  of  natural  or  man-made  targets. 
With  the  advent  of  satellite  thermal  sensor  systems  (e.g., 
the  Heat  Capacity  Mapping  Mission  and  the  proposed 
Thematic  Mapper  on  Landsat  D  as  well  as  multispectral 
aircraft  scanners),  it  is  becoming  more  important  to 
relate  underlying  scene  phenomena  to  the  remote 
sensing  observables.^  In  both  the  design  and  utilization 
of  electrooptical  sensors  it  is  important  to  be  able  to 
estimate  the  statistical  characteristics  of  the  target  and 
background  as  a  function  of  sensor  parameters.  Often 
empirical  methods  are  used  to  obtain  these  required 
data.  2*3 

An  alternative  approach  is  to  employ  a  physically 
based  or  process-oriented  model  of  the  scene.^*^  Tbis 
latter  approach  is  particularly  useful  and  often  required 
when  detailed  or  subtle  characteristics  of  a  target  need 
to  be  enhanced. 

Many  thermal  models  exist  for  different  nonvege- 
tated  targets  of  interest  and  for  planar  solid  objects. 
For  example,  Watson®  developed  a  thermal  model  for 
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predicting  the  diurnal  surface  temperature  vat .'  *'on  of 
the  ground,  and  the  University  of  Michigan^  developed 
a  model  for  the  prediction  of  time-dependent  temper¬ 
atures  and  radiance  of  planar  targets  and  backgrounds. 
However,  few  thermal  models  exist  for  plant  canopies; 
and  in  agriculture  and  forestry  applications,  vegetation 
is  the  primary  target  of  interest. 

Gates*  present^  an  energy  budget  for  a  single  plant 
leaf  isolated,  in  space,  as  did  Kimes  et  al.^  and  Wiebelt 
and  Henderson.^*  ()ther  investigators  have  modeled 
the  thermal  dynamics  of  vegetation  canopies  assuming 
a  simplistic  single  homogeneous  layer  abstraction.  For 
example,  vegetation  was  treated  as  a  single  homoge¬ 
neous  layer  with  an  associated  transmission  factor  for- 
solar  radiation  in  the  University  of  Michigan  model.'' 
Heilman  et  al.  used  thermal  scanning  data  to  measure 
crop  effective  radiant  temperatures  and  used  an  evap- 
otranspiration  (ET)  equation  to  estimate  crop  ET. 
They  assumed  that  the  sensor  was  viewing  only  the  top 
layer  of  the  crop,  and  they  ignored  the  effects  of  canopy 
geometry  on  sensor  response. 

It  is  known  that  vegetation  canopies  are  non-Lam- 
bertian  at  optical  wavelengths  primarily  due  to  canopy 
geometry.'*  Similarly,  in  the  thermal  region,  it  is  be¬ 
lieved  that  while  individual  Canopy  elements  are  iso¬ 
tropic  radiators,  the  response  from  the  canopy  may  also 
be  non-Lambertian,  because  canopy  geometry  causes 
spatial  variations  in  many  energy  flow  processes. 

A  primary  objective  of  this  paper  is  to  describe  a 
thermal  IR  exitance  model  of  a  plant  canopy,  which 
includes  a  detailed  accounting  of  vegetation  structure. 
The  model  was  applied  to  a  lodgepole  pine  canopy  over 
a  diurnal  cycle.  This  treatment  is  then  shown  to  result 
in  angular  variations  in  predicted  thermal  radiance. 
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Thermal  Vegetation  Canopy  Model  Studies 


J.  A.  SMITH,  K.  J.  RANSON,  D.  NGUYEN.  L.  BALICK, 
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An  iterative-type  thermal  model  applicable  to  forest  canopies  was  tested  with  data  from  two  diverse  forest  types.  The 
model  framewoilc  consists  of  a  system  of  steady-state  energy  budget  equations  describing  the  interactions  of  short-  and 
long-wave  radiation  within  three  horizontally  infinite  canopy  layers.  A  state-space  formulation  of  the  energy  dynamics 
within  the  canopy  is  used  which  permits  a  factorization  of  canopy  geometrical  parameters  from  cattopy  optical  and 
thermal  coefficients  as  well  as  environmental  driving  variables.  Two  sets  of  data  characterizing  a  coniferous 
(Douglas-fir)  and  deciduous  (oak~hickory)  canopy  were  collected  to  evaluate  the  thenml  model.  The  results  show  that 
the  model  approximates  measured  mean  canopy  temperatures  to  within  2°C  lor  relatively  clear  weather  conditions  and 
deviates  by  a  maximum  oi  3°C  for  very  hazy  or  foggy  conditions. 


Introduction 


Rapid  and  accurate  assessment  of 
renewable  resources  is  an  increasingly  im¬ 
portant  task  facing  remote  sensing  spe¬ 
cialists.  Mathematical  abstraction  of  en¬ 
ergy  processes  of  vegetation  canopies  is  a 
useful  technique  for  relating  sensor  re¬ 
sponse  to  environment-canopy  interac¬ 
tions.  Such  an  understanding  is  required 
in  order  to  make  timely  inferences  about 
the  condition  of  forestry  and  agricultural 
resources  from  remote  sensors.  In  the 
thermal  regime,  several  authors  have  re¬ 


ported  success  in  estimating  evapo- 
transpiration  of  crops  from  thermal  sensor 
data  (Heilman  et  al.,  1976;  Reginato 
et  al.,  1976;  and  Soer,  1980).  Several 
models  have  been  reported  that  describe 
the  energy  balance  of  vegetation  either  in 
terms  of  a  single  leaf  (Gates,  1968; 
Wiebelt  and  Henderson,  1977;  and  Kimes 
et  al.,  1978),  or  an  abstract  layered  canopy 
(Alderfer  and  Cates,  1971;  and  Deardorff, 
1978).  Few  models  have  been  described 
that  characterize  the  energy  flows  within 
vegetation  canopies  as  a  function  of 
the  canopy  geometry  (Coudriaan,  1977; 
Norman,  1979;  and  Kimes  et  al.,  1981). 
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Reflectance  of  a  vegetation  canopy  using  the  Adding  method 


K.  Cooper,  J.  A.  Smith,  and  D.  Pitts 


The  Adding  method  has  been  used  to  calculate  vegetation  canopy  bidirectional  reflectance.  Appropriate 
reflection  and  transmission  matrices  for  canopy  and  underlying  soil  layers  are  developed  in  terms  of  canopy 
geometry  and  basic  optical  properties.  The  model  is  illustrated  by  comparisons  with  field  reflectance  mea¬ 
surements  for  blue  grama  (Bouteloua  gracilis)  and  soybeans  (Calahan  9250|.  Calculations  using  the  Suits 
model  are  also  included  for  comparison. 


I.  Introduction 

The  interaction  of  electromagnetic  radiation  with 
vegetation  canopies  is  an  intriguing  problem  and  an 
understanding  of  such  interactions  has  application  to 
the  remote  sensing  of  renewable  resources.  During  the 
past  ten  years  several  investigators  have  applied  mod¬ 
eling  techniques  to  predict  vegetation  reflectance  or  as 
an  aid  in  interpreting  remotely  sensed  data.  The 
modeling  problem  generally  separates  into  two  parts. 
The  first  part  is  the  selection  of  an  appropriate  math¬ 
ematical  solution  technique  to  describe  the  radiative 
transfer  interactions  within  the  canopy.  The  second 
requirement  is  the  development  of  appropriate  ab¬ 
stractions  of  a  canopy  which  permit  the  determination 
of  electromagnetic  parameters  in  terms  of  biophysical 
attributes. 

Existing  models  which  have  been  applied  to  the  1-D 
canopy  case  include  the  Allen  et  al.,^  Suits, ^  and  Park 
and  Deering,'*  models  based  on  the  Duntley  approach, 
the  Monte  Carlo  approach  of  Smith  and  Oliver,^  the 
application  of  first-order  scattering  theory  by  Ross  and 
Nilson,®  and  the  method  of  discrete  ordinates  used  by 
Weinman  and  Guetter.^  Models  based  solely  on  geo¬ 
metric  optics  are  being  evaluated  for  selected  terrain 
materials;  for  example,  Richardson  et  al.^  Egbert,® 
Jackson,  et  al.?  and  Strahler  and  Li.‘®  Finally,  com¬ 
bined  geometric  optics  and  radiative  transfer  formu¬ 
lations  of  the  canopy  reflectance  problem  are  being  at¬ 
tempted,  must  notably  by  Verhoef  and  Bunnik"  and 
Welles  and  Norman.*^ 
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In  this  paper  we  describe  a  new  vegetation  reflectance 
model  based  on  the  Adding  method  of  van  de  Hulst.^® 
The  technique  has  been  modified  slightly  to  adapt  it  to 
the  interaction  of  shortwave  radiation  within  a  vege¬ 
tation  canopy.  The  required  scattering  matrices  for  the 
canopy  layers  are  developed  using  an  abstraction  of  the 
canopy  similar  to  that  of  the  model  of  Smith  and  Oliver; 
that  is,  reflecting  and  transmitting  leaf  facets  whose 
orientations  may  be  described  by  a  leaf-slope  distribu¬ 
tion.  Example  comparisons  of  the  model  with  field 
measurements  and  with  the  Suits  model  are  also 
given. 

II.  Adding  Method 

Consider  two  arbitrary  horizontal  layers  X,  and  X,-i 
of  a  scattering  medium  so  that  X,  lies  above  X,_i,  with 
the  following  properties:  ( 1 )  all  radiation  incident  from 
above  on  X,_i  is  either  absorbed  or  scattered  back  up¬ 
ward;  and  (2)  radiation  incident  on  X,-  may  be  re¬ 
flected,  transmitted,  or  absorbed.  Let  s,-  describe  flux 
incident  on  X,  from  above  and  define  R,  and  T,  to  be 
operators  for  layer  i  so  that  the  flux  from  s,  reflected 
back  upward  is  R,s,  and  that  transmitted  toward  layer 
I  -  1  is  T,s,.  If  Ri-\  describes  the  response  of  Xj-i  to 
incident  radiation,  the  flux  scattered  back  from  X,-i 
to  X,  is  given  by  /?,_  i(T,.s, ).  Finally,  define  two  addi¬ 
tional  operators  R;  and  T;  for  X,  which  describe  how 
layer  i  responds  to  flux  incident  from  below.  The  flux 
reflected  from  X,-- 1  and  transmitted  through  X,-  is  given 
by  Ti(R,_iT,s,),  while  that  reflected  back  toward  X,-i 
is  Ri(R,_iT,Si). 

The  scattering  of  flux  back  and  forth  between  the  two 
layers  is  indicated  schematically  in  Fig.  1  and  leads  to 
a  series  to  describe  the  total  flux  sj  escaping  X,-  up¬ 
ward; 

«;  -  iR,  +  t;«,-,t.  -I-  t;«,-,r;r,-iT,  + . .  .|»,, 

»;  -  |r,  +  t;k.-,  [  £  (R;/?.-,)"j  t.|*,. 
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Vegetation  Reflectance 
Measurements  as  a  Function  of 
Solar  Zenith  Angle 

The  wide  variability  in  diurnal  reflectance  is  caused  by  variations  in 
anisotropic  sky  irradiance,  canopy  component  geometry  and  optical 
properties,  and  type  of  reflectance  measurement. 


Introduction 

N  UNDERSTANDING  of  solar  radiation  interaction 
with  vegetation  canopies  is  necessary  for  ac¬ 
curately  interpreting  remotely  sensed  data.  Broad- 
and  narrow-band  spectral  reflectance  measure¬ 
ments  of  vegetation  canopies  are  oflen  used  to 
characterize  this  interaction.  Measured  reflec¬ 
tance,  however,  is  a  complex  function  of  canopy 


1975;  Jarvis  et  al.,  1976),  and  sen.sor  inclination 
and  azimuthal  view  angles  (Smith  and  Oliver, 
1974;  Kriehel,  1978;  Smith  et  al.,  1979). 

Understanding  variations  in  canopy  reflectance 
as  a  function  of  solar  zenith  angle  is  important  for 
Several  remote  sensing  applications.  For  example, 
such  knowledge  can  improve  multitemporal  veg¬ 
etation  classification  by  using  sun-angle  signature 


Abstract:  An  understanding  of  the  behavior  of  vegetation  canopy  reflectance 
as  a  function  of  solar  zenith  angle  is  important  to  several  remote  sensing  appli¬ 
cations.  Spectral  hemispherical-conical  reflectances  of  a  nadir  looking  sensor 
were  taken  throughout  the  day  of  a  lodgepole  pine  and  two  grass  canopies. 
Mathematical  simulations  of  both  a  spectral  hemispherical-conical  reflectance 
factor  and  a  spectral  hi-hemispherical  reflectance  were  performed  for  two 
theoretical  canopies  of  contrasting  geometric  structure.  These  results  and  re¬ 
sults  from,  literature  studies  showed  a  great  amount  of  variability  of  vegetation 
canopy  reflectances  as  a  function  of  solar  zenith  angle.  Explanations  for  this 
variability  are  discussed  and  recommendations  for  future  measurements  are 
profiosed. 


constituent  optical  properties  (Cates,  1970;  Knip- 
ling,  1970),  canopy  geometry  (Ross,  1976;  Kimes 
et  al.,  1979b)  optical  properties  of  the  ground,  at¬ 
mospheric  conditions  (KrielN!l,  1976;  Ross,  1976), 
solar  zenith  angle  (Smith  et  al.,  1974;  Kriebel, 

*  Now  with  the  Laboratory  for  Applications  of  Remote 
Sensing  (LARS),  Purdue  University,  1220  Potter  Drive, 
West  LaCsyette.  IN  47906. 


extension  techniques  (Smith  et  al.,  1975).  At 
higher  latitudes  low  sun-angles  predominate  and 
an  understanding  of  the  reflectance  changes  at  low 
sun-angles  would  be  lieneiicial.  Uiumal  ndlec- 
tance  trends  are  also  important  in  photosynthetic 
and  productivity  studies  (Kimes  et  al.,  19W). 

To  better  understand  these  relationships,  spec¬ 
tral  reflectance  measurements  were  obtained  for 
several  solar  zenith  angles  for  a  lodgepole  pine 
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Influence  of  Sky  Radiance 
Distribution  on  the  Ratio  Technique 
for  Estimating 
Bidirectional  Reflectance 

The  error  induced  in  the  estimation  of  bidirectional  reflectance 
factors  using  the  standard  ratio  technique  is  less  than  five  percent  for 
zenith  view  and  sun  angles  less  than  55  degrees. 


Introduction 

The  bidirectional  reflectance  distribution 
function,  BRDF,  of  a  surface  is  an  intrinsic 
property  of  the  material  and  is  independent  of  in¬ 
cident  irradiance  (Kriebel,  1976).  In  reality,  one 
generally  utilizes  ^e  average  of  the  BRDF  over  fi¬ 
nite  solid  angles  qf  incidence  and  exitance.  This 
average  quantity  is  termed  the  (bidirectional)  re¬ 
flectance  factor,  R,  and  is  also  defined  to  be  the 


standard  Lambertian  reference  panel  for  varying 
view  and  sun  angles. 

Several  authors  have  recently  highlighted  the 
potential  errors  induced  in  this  method  by  a  time 
varying  irradiance  field  and  have  suggested  po¬ 
tential  improvements  in  the  measurement  proce¬ 
dures  (Duggin,  1980;  Duggin,  1981;  Milton,  1981; 
Richardson,  IMl).  A  second  potential  source  of 
error  is  that  induced  by  the  difiuse  sky  radiance 
field.  The  interaction  of  the  angular  sky  irradiance 


Abstract:  The  technique  of  retioing  scene  radiance  to  the  radiance  ohtairted 
from  standard  Lambertian  reference  panels  in  order  to  estimate  bidirectional 
reflectance  factors  may  depend  on  the  angular  distribution  of  the  diffuse  <r> 
radiance  field  as  well  as  the  direct  solar  irradiance.  A  simulation  study  u?as 
performed  to  estimate  the  magnitude  of  this  effect  for  differing  dear  sky  ir- 
radiance  distributions  for  a  variety  of  vegetated  surfaces.  For  the  seven  surfaces 
and  wavelengths  atralysed,  the  error  induced  in  the  estimation  of  bidirectiorral 
reflectance  factors  using  the  standard  ratio  technique  was  less  than  5  percent 
for  xenith  view  and  sun  angles  less  than  55  degrees. 


ratio  of  the  radiant  flux  reflected  by  the  surface  to 
that  which  similarly  would  be  refitted  by  a  Lam¬ 
bertian  surface  under  the  same  viewing  and  il¬ 
lumination  geometries  (Robinson  and  Biehl, 
1979).  A  commonly  used  measurement  mefliod  to 
estimate  (bidirectional)  reflectance  factors  is  to 
ratio  the  radiance  of  a  target  to  the  radiance  of  a 
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distribution  with  target  reflectance  anisotropies 
yields  contributions  to  the  target  radiance  in  a^ 
dition  to  drat  induced  by  the  direct  solar  term. 
In  order  to  yield  a  bidiie^onal  reflectimee  iictor 
corresponding  to  the  solar  incidence  direction, 
these  contributions  of  the  sky  radiance  distribu¬ 
tion  are  implicitly  assumed  to  be  negligible.  For 
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Abstract 

Electro-optical  field  measurements  of  scene  reflectance,  together  with  appropriate  target 
characterization,  will  play  an  increasingly  Important  role  in  the  testing,  intercomparison, 
and  extension  of  scene  radiation  models.  Such  models,  in  turn,  will  serve  to  extrapolate 
the  range  of  applicability  of  field  measurements  to  a  wider  set  of  environmental  conditions 
and  sensor/source  viewing  geometries  than  could  practically  be  obtained  by  direct  measure¬ 
ment.  This  paper  presents  several  issues  which  must  be  addressed  before  such  a  symbiotic 
relationship  between  modeling  and  measurement  approaches  can  be  achieved.  The  direct  com¬ 
parison  of  model  predictions  with  field  measurements  is  seen  to  be  a  somewhat  more  subtle  task 
than  might  be  suspected. 


Introduction 


In  order  to  intelligently  design  sensor  systems  operating  in  the  optical  reflective  re¬ 
gime  to  extract  or  identify  earth  surface  terrain  features,  one  must  have  an  understanding 
of  the  basic  electro-optical  properties  of  such  features.  With  our  increased  capability  to 
measure  radiation  patterns  at  varying  spatial  and  temporal  scales,  the  limiting  factor  now 
appears  not  so  much  how  to  construct  such  devices ,  but  rather  what'  their  operating  param¬ 
eters  should  be.  In  the  remote  sensing  field,  similar  comments  are.'also  applicable  to  our 
data  processing  capabilities.  He  have  the  technological  capability  to  process  large  data 
volumes  in  a  timely  fashion,  but  still  need  to  utilize  fundamental  understanding  of  terrain 
element  behavior  in  the  design  of  appropriate  algorithms.  Smith  gave  a  recent  overview  of 
the  role  of  scene  radiation  models  in  remote  sensing  applications.^  It  is  evident  from  the 
review  that  in  the  past  20  years  several  candidate  calculational  procedures  and  media  ab¬ 
stractions  have  been  developed  to  relate  reflectance  behavior  to  intrinsic  target  charac¬ 
teristics.  In  order  to  effectively  utilize  these  candidate  tools,  the  system  designer  must 
have  a  feeling  for  their  relative  trade-offs  and  range  of  applicability  to  his  particular 
problem  areas.  To  achieve  this  goal,  the  models  must  be  subjected  to  field  seasurement 
comparisons.  The  purpose  of  this  paper  is  to  raise  some  issues  that  should  be  considered 
in  such  modeling/measurement  domparisons  and  to  point  out  some  potential  pitfalls  to  be 
avoided.  It  should  become  evident  that  field  measurement  experiments  developed  in  support 
of  modeling  activities  must  consider  additional  factors  than  those  programs  whose  primary 
mission  is  to  develop  empirical  characterizations  of  overall  terrain  element  electro-optical 
behavior. 

In  this  paper,  a  brief  overview  of  some  broad  issues  is  first  given,  followed  by  an  ex¬ 
ample  illustrating  some  of  the  issues  raised. 

Issues 


The  fundamental  intrinsic  physical  property  governing  the  reflectance  behavior  of  a  scene 
element  Is  its  bidirectional  reflectance  distribution  function,  BRDF^.  An  objective  of 
reflectance  modeling  is  typically  to  develop  relationships  between  the  BROF  and  biophysical 
or  other  meaningful  attributes  of  the  terrain  element  under  study.  In. a  strict  sense,  the 
BRDF  is  unmeasureable  and  more.comnonly,  it  is  the  intjegration  of  .this  quantity  over  finite 
solid  angles  of  incidence  and  exitance  yielding  the  reflectance  factor^  which  is  actually 
estimated  in  field  measurements.  Both  of  these  quantities  are  precisely  defined,  but  the 
fact  that,  they  vary  spectrally,  directionally,  and  that  it  is,  in  practice,  difficult  to 
develop  precise  and  dlstlct  definitions  of  target  categories,  leads  to  considerable  prag¬ 
matic  confusion  as  is  evident  in  the  literature.  In  the  natural  world,  these  quantities 
must  really  be  considered  stochastic  variables  which  means  that  one  is  always  sampling  from 
a  distribution  whose  properties  are  really  presently  poorly  understood.  Also,  since  the 
quantity  directly  measured  by  radiometers  is  the  radiance,  irtiich  is  a  convolving  of  the 
bidirectional  reflectance  distribution  function  with  the  radiance  field  and  instrument 
characteristics,  indirect  experimental  techniques,  subject  to  all  their  whimsy,  are  actual¬ 
ly  employed. 

In  this  paper  the  following  problem  is  analysed:  "What  are  some  of  the  Issues  that  must 
be  considered  when  electro-optical  field  measurements  are  obtained  in  support  of  physical 
oriented  reflectance  models  (or,  perhaps,  vice-versa)?”  Three  broad  issues  are  considered, 
to  wit : 
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Inclusion  of  a  Simple  Vegetation  Layer  in  Terrain 
Temperature  Models  for  Thermal  IR  Signature 

Prediction 

LEE  K.  BAUCK,  R.  K.  SCOGGINS,  and  L.  E.  UNK 


AteMcf-A  time-dcpeiidcnt  cMifir  budget  model  WM  dctigiicd  to  M* 
•Me  the  piedicthm  of  the  tempentore  of  tettain  teene  elementt  that 
eotitatai  a  Nmple  layet  of  vegetatioii  and  todiagnoee  the  effect  of  «cge> 
tatioR  on  lemoteljr  ien«ed  tempeiatofee.  The  model  wat  deeeloped  for 
•ee  ai  a  modnie  in  coRjunction  with  exiiting  non«egetated*leffaln  tern* 
peiatom  modelf.  Vegetation  wat  amumed  to  be  a  hoiiiontalljr  homoge- 
neont  bat  poiout  layer  partially  coveting  a  tpedfled  gtoand  larface. 
Eneigy  ba^tf  for  the  foliage  and  the  gioand  were  evalnated  eepa- 
ratcly  bat  are  interdependenL 

The  aenilthrity  of  the  vegetadon  modole  (VEGIE),  to  Ita  inpat  vari- 
dblet  when  oaed  in  tandem  with  a  ban-poond  model,  TSTM,  wae  ex¬ 
amined  and  foand  to  be  moot  ttrongly  dependent  on  the  degme  of 
fodage  coreL  The  model  reanlla  weie  verified  agdiut  meaaaienynti 
made  on  two  moderate  daya  for  a  10<m  hi^-gtam  canopy  In  Germany 
aiM  were  compared  to  reanlta  from  a  complex  vegetation  modeL  Re- 
••Its  wem  highly  aatfafaetoiy  and  ihnilar  for  both  modda.  Bodimodeia 
were  alao  applied  to  forcat  canopies  whem  aaeompdona  of  die  simple 
model  were  not  valid.  At  expected,  the  simple  model  did  not  perform 
wed,  bat  the  complex  model  wm  sadsfactoiy.  The  VEGIEnSTM  V*- 
atm  WM  then  applied  to  diagnosing  the  effects  of  varioaa  amoonts  of 
fodage  cover  over  low  emittance  tod,  with  and  withoat  leflecdon  of 
shy  thermal  IR  energy. 


I.  Introduction 

A.  Background 

Recent  analytical  and  theoretical  studies  of  the 
^themul  infrared  (IR)  emission  characteristics  of  terrain 
surfaces  have  generally  ignored  the  effects  of  vegetation  on 
thermal  IR  signatures  [1]  -[3] .  This  seems  to  be  due,  in  part, 
to  the  lack  of  a  usable  tool  for  examining  the  effects  of  vegeta¬ 
tion  in  real-world  situations.  In  truth  the  problem  is  very  com¬ 
plex  Modeb  of  vegetation  temperature  which  can  be  directly 
applied  to  remote  Knsing  problems  [4] -[6]  are  complex  and 
require  careful  spedfication  of  intracanopy  meteorological 
conditions,  canopy  structure,  and  biophysical  characteristics 
which  are  not  often  available.  Observational  and  analytical 
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studies  are  consistent  with  the  modeb  in  detailing  the  com- 
iriexity  of  the  problem  [7]  -  [  1 5  ] . 

A  major  goal  of  ongoing  work  at  the  U.S.  Army  EnginMf 
Waterways  Experiment  Station  (WES)  b  to  develop  the  ability 
to  predict  realbtic  thermal  signatures  of  natural  and  cultural 
ter^n  features  for  any  time-of-day  and  weather  conditioru. 
Woric  on  modeb  for  complex  vegetation  canopies  and  for 
planar  unvegetated  surfaces  has  pointed  to  the  need  to  fill 
the  gap  between  these  extremes.  A  model  for  predicting  the 
time-varying  temperature  for  an  unvegetated  planar  surface 
using  material  thermophysical  properties  and  meteorological 
conditions  was  developed  under  previous  WES  researdi  (16]. 
Thb  model,  the  terrain  surface  temperature  model  (TSTM),  b 
briefly  discussed  in  Part  II  of  thb  paper.  Logically,  the  next 
step  wat  to  develop  a  module,  or  submodel,  for  use  in  conjunc¬ 
tion  with  the  TSIM  to  account  for  the  dominant  effeett  of  a 
simple  byer  of  vegetation  on  thermal  IR  signatures  of  the  ter¬ 
rain  surface.  Such  a  capability  would  be  useful  even  if  it  only 
applied  to  the  simplest  of  canopies  in  nonextreme  environ¬ 
ments.  It  would  allow  the  TSTM  to  be  extended  to  ateas  of 
town,  pasture,  and  perhaps  rangelands.  The  module  developed 
in  this  context  was  named  VEGIE  and  b  described  in  thb 
paper. 

B.  Objectiva  and  Scope 

The  specific  objectives  of  the  work  presented  herein  were  to 
develop  the  capability  to  predict  the  temperature  of  terrain 
scene  elements  that  contain  a  simple  layer  of  vegetation  and  to 
dbgnose  the  effect  of  vegetation  on  remotely  sensed  tempera¬ 
tures  of  terrain  elements.  VEGIE  was  designed  for  immedbte 
application  to  thermal  IR  signature  prediction  and  analytical 
studies  while  more  complete  and  theoretical  treatments  of  the 
probbm  were  under  way.  It  was  also  required  that  VEGIE 
maintain  the  operational  flexibility  and  sim^icity  of  the  planar 
surface  nsodel,  TSTM.  Some  of  the  more  important  of  the 
TSTM  eharMterbtics  are: 

0  time  dependence  and  fast  response  to  environmental 
changes; 

2)  air  temperature  considered  a  state  variable; 

3)  materlab  treated  as  horizontally  and  vertically  homo¬ 
geneous  layers; 


0I96-2892/8I/0700-0I43$00.7S  0 1981  IEEE 


AlO 


Radiative  transfer  modei  for  heterogeneous  3-D  scenes 

D.  S.  Kimes  and  J.  A.  Kirchner 


A  general  mathematical  framework  simulating  processes  in  heterogeneous  3-D  scenes  is  presented.  Spe- 
cincally,  a  model  was  designed  and  coded  for  application  to  radiative  transfers  in  vegetative  scenes.  The 
model  is  unique  in  that  it  predicts  (1)  the  directional  spectral  reflectance  factors  as  a  function  of  the  sensor’s 
azimuth  and  zenith  angles  and  the  sensor’s  p<»ition  above  the  canopy,  (2)  the  spectral  alisorption  as  a  func¬ 
tion  of  location  within  the  scene,  and  (3)  the  directionBl  spectral  radiance  as  a  function  of  the  sensor’s  loca¬ 
tion  within  the  scene.  The  model  was  shown  to  follow  known  physical  principles  of  radiative  transfer.  Ini¬ 
tial  veriflcation  of  the  model  as  applied  to  a  soybean  row  crop  showed  that  the  simulated  directional  reflec¬ 
tance  data  corresponded  relatively  well  in  gross  trends  to  the  measured  data.  However,  the  model  can  be 
greatly  improved  by  incorporating  more  sopbisticated  and  realistic  anisotropic  scattering  algorithms. 


I.  Introductiofi 

Understanding  the  nature  of  the  interactions  between 
electromagnetic  radiation  and  terrestrial  scenes  is 
necessary  to  utilize  remotely  sensed  data  effectively. 
Mathematical  models  which  simulate  these  radiant 
interactioits  are  valuable  quantitative  tools  for  relating 
the  physical  features  of  the  terrestrial  scene  and  radiant 
sources  to  the  radiant  response  of  the  scene.  In  some 
cases  these  models  can  be  used  to  infer  underlying  scene 
attributes  using  remotely  sensed  data.  Presently,  most 
models  are  one-dimensional  in  that  they  treat  one  or 
more  layers  of  spatially  homogeneous  material  that 
extends  infinitely  in  the  horizontal  direction. 

There  is  a  need  for  a  general  modeling  framework  for 
complex  scenes  which  are  spatially  heterogeneous  in 
three  dimensions.  A  number  of  important  scenes  are 
extremely  variable  spatially  and  cannot  be  treated  by 
1-D  models.  For  example,  forest  canopies  often  have 
large  natural  openings  between  the  tree  crowns;  in¬ 
complete  row  crop  canopies  show  an  extreme  spatial 
variation  between  rows;  military  detection  and  cam¬ 
ouflage  problems  within  clumps  of  vegetation  are  3-D 
problems;  and  urban  scenes  are  also  complex  3-D 
scenes. 

The  difficulty  in  modeling  any  complex  3-D  structure 
is  accounting  for  the  many  phenomena  that  occur  with 
variable  intensities  throughout  the  scene.  The  fol¬ 
lowing  proposed  3-D  model  may  provide  a  general  ac¬ 
counting  framework  for  many  3-D  modeling  problems 
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where  processes  vary  in  intensity  and  direction  of 
movement  as  a  function  of  spatial  location.  Specifi¬ 
cally,  the  model  described  in  this  study  was  designed 
and  coded  for  application  to  radiative  transfers  in  veg¬ 
etative  scenes.  However,  the  basic  principles  of  the 
model  may  provide  a  general  framework  for  other  3-D 
scenes  and  other  processes. 

Several  mathematical  models  have  been  developed 
for  vegetation  canopies  in  the  optical  wavelengths. 
Most  of  these  models  treat  homogeneous  canopies 
which  vary  structurally  only  with  height  (1-D  models) 
and  are  reviewed  by  Cooper  et  al.  ‘  Exceptions  include 
geometric  optic  models  for  row  cit^  developed  by  both 
Richardson  et  al.^  and  Jackson  et  al.^  and  for  forest 
canopies  by  Strahler  and  Li.*  In  addition,  several 
models  combine  both  radiative  transfer  and  geometric 
optic  techniques  as  are  the  extended  versions  of  the  1-D 
Suits  model^  for  row  crops  developed  by  Verhoef  and 
Bunnik,*  the  more  general  extension  of  the  Suits  1-D 
model  for  row  crops  by  Suits,''  and  the  3-D  reflectance 
model  based  on  an  array  of  ellipsoids  by  Welles  and 
Norman.*  Several  other  models  with  more  simplistic 
assumptions  and/or  empirical  characterizations  are 
reviewed  by  Smith  and  Ranson.* 

In  this  paper  we  present  the  concepts  of  a  general 
mathematical  framework  for  simulating  processes  in 
heterogeneous  3-D  scenes.  Speciflcally,  in  this  study 
the  m^el  is  applied  to  radiant  transfers  (reflective 
wavelengths— visible  and  near  IR)  in  vegetation  cano¬ 
pies  wiU)  both  homogeneous  structure  and  row  crop 
structure.  The  model  simulates  the  multiple  scattering 
which  occurs  within  the  scene  given  the  global  distri¬ 
bution  of  solar  irradiance  conditions  and  the  optical 
properties  and  spatial  structure  of  the  scene  compo¬ 
nents.  Tlie  model  predicts  (1)  Uie  directional  spectral 
reflectance  factors  as  a  function  of  the  sensor's  azimuUi 
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